Merge branch 'par_rendering' into 'master'

Draft: Parallel rendering

See merge request vtk/vtk-m!2966

README

@@ -0,0 +1,142 @@
+vtkh:
+tests/vtkh/vtk-h_render.cpp
+
+TEST(vtkh_render, vtkh_bg_color)
+
+
+make data
+bounds
+camera setup.
+vtkh::MakeRender(...)
+   vtkh::Render: info needed to create a single image. N domains = N canvases
+     camera, bounds, FG/BG colors
+
+vtkh::RayTracer tracer;
+tracer.SetInput(...);
+
+vtkh::Scene scene;
+scene.AddRender(render);
+scene.AddRenderer(&tracer);
+scene.Render();
+
+vtkh::RayTracer : public vtkh::Renderer
+ --it returns a vtkm::rendering::CanvasRayTracer(width,height)
+ DoExecute() does the rendering
+ PostExecute() does the compositing.
+
+
+vtkh::Scene
+  list of vtkh::Renderer (raytracer, ...)
+  vector of vtkh::Render (camera, data, fg/bg)
+  Scene::Render():
+    for each renderer:
+      renderer->Update() does the rendering.
+    (do opaque, then volume)
+ RenderWorldAnnotations()
+   Render::RenderWorldAnnotations()
+     if rank != 0 then RETURN
+     Annotator annotator(canvas, camera, m_scene_bounds)
+     annotator.RenderWorldAnnotations()
+ RenderScreenAnnotations()
+
+vtkh::Image. pixels, depth, composite order, ...
+
+
+vtkh::Renderer : public Filter
+  DoExecute:
+     for each ds in m_input
+       m_mapper->RenderCells(ds)
+     this->Composite(...)
+
+
+vtkh::Compositor
+  Composite()
+    calls RadixKCompositor, etc.
+
+
+
+vtkh::Render
+  camera
+  image name
+  width/height/bounds
+  fg/bg colors
+
+vtkh::RayTracer
+  SetInput: dataset
+  SetField: scalar
+
+vtkh::Scene
+  AddRender: vtkh::Render
+  AddRenderer: vtkh::RayTracer
+
+
+
+vtkh::Scene.Render();
+
+
+
+=======================================================================================
+vtkm:
+vtkm::rendering::Scene
+  AddActor: dataset
+  Render(mapper, canvas, camera)
+
+vtkm::rendering::View
+  camera
+  canvas
+  mapper
+  scene
+  annotations
+  Paint() does the rendering
+
+vtkm::rendering::Canvas
+ Fg/Bg
+ color/depth buffers
+ model/view mtx
+
+
+============================================
+VTK-m: simpler example... ?
+examples/demo/Demo.cxx
+
+vtkm::rendering::Camera (look at, up, clip, etc).
+vtkm::rendering::Actor (dataset)
+vtkm::rendering::CanvasRayTracer (x,y)
+
+vtkm::rendering::Scene (actor)
+vtkm::rendering::View3D view (scene, mapper, canvas, camera, bg)
+
+view.Paint() (renders the image).
+
+
+
+
+
+=========================================================================================
+Comparing classes, etc.
+Camera: the same. vtkm::rendering::Camera
+
+Scene:
+vtkm: vector of vtkm::rendering::Actor
+vtkh: list/vector of vtkh::Renderer, batchsize, has_volume.
+
+vtkh::Render aprox equal to vtkm::rendering::View
+
+
+
+
+
+
+==============================================================================
+vtkh volume render unstructured:
+
+
+scene.AddRender()
+scene.AddRenderer(vtkh::VolumeRenderer);
+scene.Render()
+  // pass 2 for volume
+  renderer->SetComponosite(true)
+  renderer->SetRenders(current_batch);
+  renderer->Update() //does the rendering...
+     DoExecute()
+       vtkm::rendering::Mapper::RenderCells()

examples/CMakeLists.txt

@@ -52,6 +52,8 @@ if(VTKm_ENABLE_EXAMPLES)
   add_subdirectory(temporal_advection)
   add_subdirectory(tetrahedra)
   add_subdirectory(smoke_test)
+
+  add_subdirectory(frontier_benchmark)
 endif()
 
 if (VTKm_ENABLE_TESTING)

examples/demo/Demo.cxx

@@ -8,6 +8,8 @@
 //  PURPOSE.  See the above copyright notice for more information.
 //============================================================================
 
+#include <vtkm/io/VTKDataSetWriter.h>
+
 #include <vtkm/cont/Initialize.h>
 #include <vtkm/source/Tangle.h>
 
@@ -39,6 +41,9 @@ int main(int argc, char* argv[])
   vtkm::cont::DataSet tangleData = tangle.Execute();
   std::string fieldName = "tangle";
 
+  vtkm::io::VTKDataSetWriter writer("tangle.vtk");
+  writer.WriteDataSet(tangleData);
+
   // Set up a camera for rendering the input data
   vtkm::rendering::Camera camera;
   camera.SetLookAt(vtkm::Vec3f_32(0.5, 0.5, 0.5));
@@ -57,7 +62,7 @@ int main(int argc, char* argv[])
   vtkm::rendering::Scene scene;
   scene.AddActor(actor);
   // 2048x2048 pixels in the canvas:
-  CanvasRayTracer canvas(2048, 2048);
+  CanvasRayTracer canvas(512, 512);
   // Create a view and use it to render the input data using OS Mesa
 
   vtkm::rendering::View3D view(scene, MapperVolume(), canvas, camera, bg);

examples/frontier_benchmark/CMakeLists.txt

@@ -0,0 +1,19 @@
+##============================================================================
+##  Copyright (c) Kitware, Inc.
+##  All rights reserved.
+##  See LICENSE.txt for details.
+##
+##  This software is distributed WITHOUT ANY WARRANTY; without even
+##  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+##  PURPOSE.  See the above copyright notice for more information.
+##============================================================================
+cmake_minimum_required(VERSION 3.12...3.15 FATAL_ERROR)
+project(VTKmFrontierBenchmark CXX)
+
+#Find the VTK-m package
+find_package(VTKm REQUIRED QUIET)
+
+if(TARGET vtkm::rendering AND TARGET vtkm::filter_contour AND TARGET vtkm::source)
+  add_executable(FrontierBenchmark FrontierBenchmark.cxx)
+  target_link_libraries(FrontierBenchmark PRIVATE vtkm::rendering vtkm::filter_contour vtkm::source)
+endif()

examples/frontier_benchmark/FrontierBenchmark.cxx

@@ -0,0 +1,457 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#include <vtkm/Bounds.h>
+#include <vtkm/cont/Timer.h>
+#include <vtkm/filter/contour/Contour.h>
+#include <vtkm/rendering/Actor.h>
+#include <vtkm/rendering/CanvasRayTracer.h>
+#include <vtkm/rendering/MapperRayTracer.h>
+#include <vtkm/rendering/Scene.h>
+#include <vtkm/rendering/View3D.h>
+#include <vtkm/rendering/compositing/Image.h>
+#include <vtkm/rendering/testing/RenderTest.h>
+#include <vtkm/source/PerlinNoise.h>
+
+#include <iomanip>
+#include <iostream>
+#include <random>
+#include <sstream>
+
+const static std::string PERLIN_MODE_GROW = "grow";
+const static std::string PERLIN_MODE_SUBDIVIDE = "subdivide";
+const static std::string CAMERA_MODE_STATIC = "static";
+const static std::string CAMERA_MODE_ORBIT = "orbit";
+
+struct BenchmarkOptions
+{
+  BenchmarkOptions(int argc, char** argv) { this->Parse(argc, argv); }
+
+  void Parse(int argc, char** argv)
+  {
+    for (int i = 1; i < argc; ++i)
+    {
+      auto arg = std::string(argv[i]);
+      auto seperatorPos = arg.find('=');
+      if (seperatorPos == std::string::npos)
+      {
+        seperatorPos = arg.size();
+      }
+
+      std::string key = arg.substr(0, seperatorPos);
+      std::string val = "";
+      if (seperatorPos != std::string::npos && seperatorPos + 1 < arg.size())
+      {
+        val = arg.substr(seperatorPos + 1);
+      }
+
+      if (key == "--perlin-dims")
+      {
+        int start = 0;
+        for (int d = 0; d < 3; ++d)
+        {
+          auto end = val.find(',', start);
+          auto dim = val.substr(start, end - start);
+          this->PerlinDimensions[d] = std::stoi(dim);
+          start = end + 1;
+        }
+      }
+      else if (key == "--perlin-seed")
+      {
+        this->PerlinSeed = std::stoi(val);
+      }
+      else if (key == "--perlin-mode")
+      {
+        this->PerlinMode = val;
+      }
+      else if (key == "--perlin-scale")
+      {
+        this->PerlinScale = std::stof(val);
+      }
+      else if (key == "--width")
+      {
+        this->CanvasWidth = std::stoi(val);
+      }
+      else if (key == "--height")
+      {
+        this->CanvasHeight = std::stoi(val);
+      }
+      else if (key == "--iters")
+      {
+        this->NumIterations = std::stoi(val);
+      }
+      else if (key == "--timing-file")
+      {
+        this->TimingFileName = val;
+      }
+      else if (key == "--camera-mode")
+      {
+        this->CameraMode = val;
+      }
+      else if (key == "--image-format")
+      {
+        this->ImageFormat = val;
+      }
+      else if (key == "--show-args")
+      {
+        this->ShowArgs = true;
+      }
+    }
+  }
+  vtkm::Id3 PerlinDimensions = vtkm::Id3(1024, 1024, 1024);
+  vtkm::IdComponent PerlinSeed = 1;
+  std::string PerlinMode = PERLIN_MODE_GROW;
+  vtkm::Float32 PerlinScale = 3.0f;
+  vtkm::Id CanvasWidth = 1920;
+  vtkm::Id CanvasHeight = 1080;
+  vtkm::Id NumIterations = 10;
+  std::string TimingFileName = "timing.csv";
+  std::string ImageFormat = "png";
+  std::string CameraMode = CAMERA_MODE_STATIC;
+  bool ShowArgs = false;
+};
+
+struct MpiTopology
+{
+  MpiTopology(int rank, int size)
+    : Rank(rank)
+    , Size(size)
+  {
+  }
+
+  int Rank;
+  int Size;
+
+  int XRank;
+  int YRank;
+  int ZRank;
+  int XSize;
+  int YSize;
+  int ZSize;
+
+  void SetShapeToCube()
+  {
+    int sizeCbrt = static_cast<int>(std::cbrt(static_cast<float>(this->Size)));
+    this->ZSize = sizeCbrt;
+    this->YSize = sizeCbrt;
+    this->XSize = this->Size / (this->YSize * this->ZSize);
+
+    this->XRank = this->Rank / (this->YSize * this->ZSize);
+    this->YRank = (this->Rank - (this->XRank * this->YSize * this->ZSize)) / this->ZSize;
+    this->ZRank =
+      this->Rank - (this->XRank * this->YSize * this->ZSize) - (this->YRank * this->ZSize);
+  }
+};
+
+struct IterationTimes
+{
+  vtkm::Float64 RenderTime = -1.0f;
+  vtkm::Float64 CompositeTime = -1.0f;
+  vtkm::Float64 TotalTime = -1.0f;
+};
+
+std::string GetImageName(const std::string& prefix,
+                         const BenchmarkOptions& options,
+                         const MpiTopology& mpiTopology)
+{
+  std::stringstream ss;
+  ss << options.CameraMode << "/" << prefix << "_" << options.PerlinMode << "_"
+     << options.CanvasWidth << "x" << options.CanvasHeight << "_" << mpiTopology.Size << "."
+     << options.ImageFormat;
+  return ss.str();
+}
+
+std::string GetFrameName(const std::string& prefix,
+                         int frameNumber,
+                         const BenchmarkOptions& options,
+                         const MpiTopology& mpiTopology)
+{
+  std::stringstream ss;
+  ss << options.CameraMode << "/" << prefix << "_" << options.PerlinMode << "_" << mpiTopology.Size
+     << "_" << options.CameraMode << "_frame_" << std::setw(4) << std::setfill('0') << frameNumber
+     << "." << options.ImageFormat;
+  return ss.str();
+}
+
+void CollectTimeStats(vtkm::Float64 time,
+                      vtkm::Float64& minTime,
+                      vtkm::Float64& maxTime,
+                      vtkm::Float64& avgTime)
+{
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+  std::vector<vtkm::Float64> allTimes(comm.size(), -1.0f);
+  if (comm.rank() == 0)
+    vtkmdiy::mpi::gather(comm, time, allTimes, 0);
+  else
+    vtkmdiy::mpi::gather(comm, time, 0);
+  if (comm.rank() == 0)
+  {
+    auto minMaxTime = std::minmax_element(allTimes.begin(), allTimes.end());
+    minTime = *(minMaxTime.first);
+    maxTime = *(minMaxTime.second);
+    avgTime = std::accumulate(allTimes.begin(), allTimes.end(), 0.0) /
+      static_cast<vtkm::Float64>(allTimes.size());
+  }
+}
+
+void SaveTimeStats(const std::vector<IterationTimes>& stats,
+                   const BenchmarkOptions& options,
+                   const MpiTopology& vtkmNotUsed(mpiTopology))
+{
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+  if (comm.rank() != 0)
+    return;
+
+  std::vector<vtkm::Float64> renderTimes;
+  std::vector<vtkm::Float64> compositeTimes;
+  std::vector<vtkm::Float64> totalTimes;
+
+  for (const auto& stat : stats)
+  {
+    renderTimes.push_back(stat.RenderTime);
+    compositeTimes.push_back(stat.CompositeTime);
+    totalTimes.push_back(stat.TotalTime);
+  }
+
+  auto CalculateTimeStats = [](const std::vector<vtkm::Float64>& times,
+                               vtkm::Float64& minTime,
+                               vtkm::Float64& maxTime,
+                               vtkm::Float64& avgTime) {
+    auto minMaxTime = std::minmax_element(times.begin(), times.end());
+    minTime = *(minMaxTime.first);
+    maxTime = *(minMaxTime.second);
+    avgTime =
+      std::accumulate(times.begin(), times.end(), 0.0) / static_cast<vtkm::Float64>(times.size());
+  };
+
+  vtkm::Float64 minRTime, maxRTime, avgRTime;
+  vtkm::Float64 minCTime, maxCTime, avgCTime;
+  vtkm::Float64 minTTime, maxTTime, avgTTime;
+  CalculateTimeStats(renderTimes, minRTime, maxRTime, avgRTime);
+  CalculateTimeStats(compositeTimes, minCTime, maxCTime, avgCTime);
+  CalculateTimeStats(totalTimes, minTTime, maxTTime, avgTTime);
+
+  auto ToHumanReadableMS = [](vtkm::Float64 time) { return static_cast<int>(time * 1000); };
+  std::ofstream file;
+  file.open(options.TimingFileName, std::ios::app);
+  bool hasHeader = file.tellp() != 0;
+  if (!hasHeader)
+  {
+    file << "World Size,Canvas Size,Min Render Time,Max Render Time,Avg Render Time,"
+            "Min Composite Time,Max Composite Time,Avg Composite Time,"
+            "Min Total Time,Max Total Time,Avg Total Time"
+         << std::endl;
+  }
+  file << comm.size() << "," << options.CanvasWidth << "x" << options.CanvasHeight << ",";
+  file << ToHumanReadableMS(minRTime) << "," << ToHumanReadableMS(maxRTime) << ","
+       << ToHumanReadableMS(avgRTime) << ",";
+  file << ToHumanReadableMS(minCTime) << "," << ToHumanReadableMS(maxCTime) << ","
+       << ToHumanReadableMS(avgCTime) << ",";
+  file << ToHumanReadableMS(minTTime) << "," << ToHumanReadableMS(maxTTime) << ","
+       << ToHumanReadableMS(avgTTime) << std::endl;
+  file.close();
+}
+
+void GenerateDataSet(const BenchmarkOptions& options,
+                     const MpiTopology& mpiTopology,
+                     vtkm::cont::DataSet& dataSet,
+                     std::string& fieldName,
+                     vtkm::Range& globalFiendRange,
+                     vtkm::Bounds& globalBounds)
+{
+  fieldName = "perlinnoise";
+  vtkm::source::PerlinNoise perlin;
+  perlin.SetCellDimensions(options.PerlinDimensions);
+  perlin.SetSeed(options.PerlinSeed);
+
+  // Perlin Noise does not generate "interesting" surfaces when the scale is too small, .i.e,
+  // X, Y, Z values are too close to each other. Hence we scale the perlin noise by a factor
+  if (options.PerlinMode == PERLIN_MODE_GROW)
+  {
+    vtkm::Vec3f origin = vtkm::Vec3f_32{ static_cast<vtkm::Float32>(mpiTopology.XRank),
+                                         static_cast<vtkm::Float32>(mpiTopology.YRank),
+                                         static_cast<vtkm::Float32>(mpiTopology.ZRank) } *
+      options.PerlinScale;
+    vtkm::Vec3f maxExtent = origin + vtkm::Vec3f{ options.PerlinScale };
+    perlin.SetOrigin(origin);
+    perlin.SetMaxExtent(maxExtent);
+    dataSet = perlin.Execute();
+  }
+  else if (options.PerlinMode == PERLIN_MODE_SUBDIVIDE)
+  {
+    vtkm::Vec3f_32 blockSpacing = vtkm::Vec3f_32{ options.PerlinScale } *
+      vtkm::Vec3f_32{ 1.0f / static_cast<vtkm::Float32>(mpiTopology.XSize),
+                      1.0f / static_cast<vtkm::Float32>(mpiTopology.YSize),
+                      1.0f / static_cast<vtkm::Float32>(mpiTopology.ZSize) };
+    vtkm::Vec3f_32 origin = vtkm::Vec3f_32{ static_cast<vtkm::Float32>(mpiTopology.XRank),
+                                            static_cast<vtkm::Float32>(mpiTopology.YRank),
+                                            static_cast<vtkm::Float32>(mpiTopology.ZRank) } *
+      blockSpacing;
+    vtkm::Vec3f_32 maxExtent = origin + blockSpacing;
+
+    perlin.SetOrigin(origin);
+    perlin.SetMaxExtent(maxExtent);
+    dataSet = perlin.Execute();
+  }
+
+  std::vector<vtkm::Float64> isoValues{ 0.4f, 0.75f };
+  vtkm::filter::contour::Contour contour;
+  contour.SetIsoValues(isoValues);
+  contour.SetActiveField(fieldName);
+  dataSet = contour.Execute(dataSet);
+
+  globalFiendRange = { 0.0f, 1.0f };
+  if (options.PerlinMode == PERLIN_MODE_GROW)
+  {
+    globalBounds = { vtkm::Vec3f{ 0.0f, 0.0f, 0.0f },
+                     vtkm::Vec3f{ static_cast<vtkm::FloatDefault>(mpiTopology.XSize),
+                                  static_cast<vtkm::FloatDefault>(mpiTopology.YSize),
+                                  static_cast<vtkm::FloatDefault>(mpiTopology.ZSize) } *
+                       options.PerlinScale };
+  }
+  else if (options.PerlinMode == PERLIN_MODE_SUBDIVIDE)
+  {
+    globalBounds = { vtkm::Vec3f{ 0.0f }, vtkm::Vec3f{ options.PerlinScale } };
+  }
+
+  // Add a small epsilon to the bounds to prevent the world annotations from being clipped
+  /*
+  vtkm::Float64 boundsEps = 0.0f;
+  if (options.CameraMode == CAMERA_MODE_ORBIT)
+  {
+    boundsEps = 0.2f;
+  }
+  globalBounds.Include(globalBounds.MinCorner() -
+                       vtkm::Vec3f_64{ boundsEps, boundsEps, boundsEps });
+  globalBounds.Include(globalBounds.MaxCorner() +
+                       vtkm::Vec3f_64{ boundsEps, boundsEps, boundsEps });
+  */
+}
+
+void RunBenchmark(const BenchmarkOptions& options)
+{
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+  MpiTopology mpiTopology(comm.rank(), comm.size());
+  mpiTopology.SetShapeToCube();
+
+  // Generate DataSet
+  vtkm::cont::DataSet dataSet;
+  std::string fieldName;
+  vtkm::Range globalFiendRange;
+  vtkm::Bounds globalBounds;
+  GenerateDataSet(options, mpiTopology, dataSet, fieldName, globalFiendRange, globalBounds);
+
+  vtkm::rendering::Scene scene;
+  vtkm::cont::ColorTable colorTable("inferno");
+  vtkm::rendering::Actor actor(
+    dataSet.GetCellSet(), dataSet.GetCoordinateSystem(), dataSet.GetField(fieldName), colorTable);
+  actor.SetScalarRange(globalFiendRange);
+  scene.AddActor(actor);
+
+  vtkm::rendering::Camera camera;
+  camera.ResetToBounds(globalBounds);
+
+  vtkm::rendering::CanvasRayTracer canvas(options.CanvasWidth, options.CanvasHeight);
+
+  if (options.CameraMode == CAMERA_MODE_STATIC)
+  {
+    camera.Azimuth(10.0f);
+    camera.Elevation(20.0f);
+    std::vector<IterationTimes> benchmarkTimes;
+    for (int iter = 0; iter < options.NumIterations; iter++)
+    {
+      vtkm::rendering::Color bg(0.2f, 0.2f, 0.2f, 1.0f);
+      vtkm::rendering::View3D view(scene, vtkm::rendering::MapperRayTracer(), canvas, camera, bg);
+      view.Paint();
+      if (comm.rank() == 0)
+      {
+        benchmarkTimes.push_back(
+          IterationTimes{ .RenderTime = view.GetTimes()[vtkm::rendering::RENDER_TIME_KEY],
+                          .CompositeTime = view.GetTimes()[vtkm::rendering::COMPOSITE_TIME_KEY],
+                          .TotalTime = view.GetTimes()[vtkm::rendering::TOTAL_TIME_KEY] });
+      }
+    }
+    SaveTimeStats(benchmarkTimes, options, mpiTopology);
+    if (mpiTopology.Rank == 0)
+    {
+      canvas.SaveAs(GetImageName("perlin_static", options, mpiTopology));
+    }
+  }
+  else if (options.CameraMode == CAMERA_MODE_ORBIT)
+  {
+    std::random_device dev;
+    std::mt19937 rng(dev());
+    std::uniform_real_distribution<vtkm::Float64> dist(0.0, 1.0);
+    vtkm::Float64 dirX = -1.0f;
+    for (int iter = 0; iter < options.NumIterations; iter++)
+    {
+      vtkm::rendering::Color bg(0.2f, 0.2f, 0.2f, 1.0f);
+      vtkm::rendering::View3D view(scene, vtkm::rendering::MapperRayTracer(), canvas, camera, bg);
+      view.Paint();
+
+      if (mpiTopology.Rank == 0)
+      {
+        canvas.SaveAs(GetFrameName("perlin_movie", iter, options, mpiTopology));
+      }
+
+      vtkm::Float64 speedX = 0.01f * dirX;
+      vtkm::Float64 speedY = 0.0f;
+      camera.TrackballRotate(0.0, 0.0, speedX, speedY);
+
+      if (mpiTopology.Rank == 0 && iter > 0 && (iter + 1) % 10 == 0)
+      {
+        std::cerr << "Frame " << (iter + 1) << " of " << options.NumIterations << " done"
+                  << std::endl;
+      }
+    }
+  }
+}
+
+int main(int argc, char* argv[])
+{
+  // Initialize MPI
+  std::unique_ptr<vtkmdiy::mpi::environment> mpiEnv = nullptr;
+  if (!vtkmdiy::mpi::environment::initialized())
+  {
+    mpiEnv.reset(new vtkmdiy::mpi::environment(argc, argv));
+  }
+
+  // Initialize VTK-m
+  vtkm::cont::Initialize(argc, argv, vtkm::cont::InitializeOptions::None);
+
+  BenchmarkOptions options(argc, argv);
+  if (options.ShowArgs)
+  {
+    std::cerr << std::boolalpha;
+    std::cerr << argv[0] << ":" << std::endl;
+    std::cerr << "\tPerlin Dimensions: " << options.PerlinDimensions << std::endl;
+    std::cerr << "\tPerlin Seed: " << options.PerlinSeed << std::endl;
+    std::cerr << "\tCanvas Width: " << options.CanvasWidth << std::endl;
+    std::cerr << "\tCanvas Height: " << options.CanvasHeight << std::endl;
+    std::cerr << "\tNum Iterations: " << options.NumIterations << std::endl;
+    std::cerr << "\tTiming File: " << options.TimingFileName << std::endl;
+    std::cerr << "\tCamera Mode: " << options.CameraMode << std::endl;
+    std::cerr << "\tShow Args: " << options.ShowArgs << std::endl;
+    std::cerr << std::noboolalpha;
+  }
+
+  RunBenchmark(options);
+  return 0;
+}

vtkm/rendering/CMakeLists.txt

@@ -55,6 +55,10 @@ set(headers
   View3D.h
   Wireframer.h
   WorldAnnotator.h
+
+  compositing/Compositor.h
+  compositing/Image.h
+  compositing/PNGEncoder.h
   )
 
 set(sources
@@ -86,6 +90,15 @@ set(sources
   raytracing/Logger.cxx
   raytracing/MeshConnectivityContainers.cxx
   raytracing/TriangleExtractor.cxx
+
+  compositing/Compositor.cxx
+  compositing/DirectSendCompositor.cxx
+  compositing/Image.cxx
+  compositing/PartialCompositor.cxx
+  compositing/PNGEncoder.cxx
+  compositing/RadixKCompositor.cxx
+  compositing/PayloadCompositor.cxx
+  compositing/PayloadImage.cxx
   )
 
 # This list of sources has code that uses devices and so might need to be
@@ -147,6 +160,10 @@ if(UNIX AND NOT APPLE)
   target_link_libraries(vtkm_rendering PRIVATE rt)
 endif()
 
+if (VTKm_ENABLE_MPI)
+  target_link_libraries(vtkm_rendering PUBLIC MPI::MPI_CXX)
+endif()
+
 #-----------------------------------------------------------------------------
 add_subdirectory(internal)
 add_subdirectory(raytracing)

vtkm/rendering/Canvas.cxx

@@ -12,6 +12,7 @@
 
 #include <vtkm/cont/ArrayHandleCounting.h>
 #include <vtkm/cont/DataSetBuilderUniform.h>
+#include <vtkm/cont/EnvironmentTracker.h>
 #include <vtkm/cont/TryExecute.h>
 #include <vtkm/io/DecodePNG.h>
 #include <vtkm/io/EncodePNG.h>
@@ -195,6 +196,27 @@ struct DrawColorBar : public vtkm::worklet::WorkletMapField
   bool Horizontal;
 }; // struct DrawColorBar
 
+struct CopyFromBuffers : public vtkm::worklet::WorkletMapField
+{
+  using ControlSignature = void(FieldInOut, FieldInOut, WholeArrayIn, WholeArrayIn);
+  using ExecutionSignature = void(_1, _2, _3, _4, WorkIndex);
+
+  template <typename ColorPortalType, typename DepthPortalType>
+  VTKM_EXEC void operator()(vtkm::Vec4f_32& color,
+                            vtkm::Float32& depth,
+                            const ColorPortalType& colorBuffer,
+                            const DepthPortalType& depthBuffer,
+                            const vtkm::Id& index) const
+  {
+    vtkm::Id colorOffset = index * 4;
+    for (vtkm::IdComponent i = 0; i < 4; ++i)
+    {
+      color[i] = static_cast<vtkm::Float32>(colorBuffer.Get(colorOffset + i) / 255.0f);
+    }
+    depth = static_cast<vtkm::Float32>(depthBuffer.Get(index));
+  }
+}; // struct CopyFromBuffers
+
 } // namespace internal
 
 struct Canvas::CanvasInternals
@@ -615,6 +637,13 @@ void Canvas::SetViewToScreenSpace(const vtkm::rendering::Camera& vtkmNotUsed(cam
 
 void Canvas::SaveAs(const std::string& fileName) const
 {
+  //Only rank 0 has the composited image.
+#ifdef VTKM_ENABLE_MPI
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+  if (comm.rank() != 0)
+    return;
+#endif
+
   this->RefreshColorBuffer();
   ColorBufferType::ReadPortalType colorPortal = GetColorBuffer().ReadPortal();
   vtkm::Id width = GetWidth();
@@ -661,5 +690,13 @@ vtkm::rendering::WorldAnnotator* Canvas::CreateWorldAnnotator() const
 {
   return new vtkm::rendering::WorldAnnotator(this);
 }
+
+void Canvas::CopyFrom(const vtkm::cont::ArrayHandle<unsigned char>& colorBuffer,
+                      const vtkm::cont::ArrayHandle<vtkm::Float32>& depthBuffer)
+{
+  vtkm::worklet::DispatcherMapField<internal::CopyFromBuffers> dispatcher(
+    internal::CopyFromBuffers{});
+  dispatcher.Invoke(this->GetColorBuffer(), this->GetDepthBuffer(), colorBuffer, depthBuffer);
+}
 }
 } // vtkm::rendering

vtkm/rendering/Canvas.h

@@ -225,6 +225,10 @@ public:
   VTKM_CONT
   void EndTextRenderingBatch() const;
 
+  VTKM_CONT
+  void CopyFrom(const vtkm::cont::ArrayHandle<unsigned char>& colorBuffer,
+                const vtkm::cont::ArrayHandle<vtkm::Float32>& depthBuffer);
+
   friend class AxisAnnotation2D;
   friend class ColorBarAnnotation;
   friend class ColorLegendAnnotation;

vtkm/rendering/MapperRayTracer.cxx

@@ -153,6 +153,11 @@ void MapperRayTracer::SetShadingOn(bool on)
   this->Internals->Shade = on;
 }
 
+void MapperRayTracer::SetLightPosition(const vtkm::Vec3f_32& lightPosition)
+{
+  // this->Internals->Tracer.SetLightPosition(lightPosition);
+}
+
 vtkm::rendering::Mapper* MapperRayTracer::NewCopy() const
 {
   return new vtkm::rendering::MapperRayTracer(*this);

vtkm/rendering/MapperRayTracer.h

@@ -39,6 +39,8 @@ public:
   vtkm::rendering::Mapper* NewCopy() const override;
   void SetShadingOn(bool on);
 
+  void SetLightPosition(const vtkm::Vec3f_32& lightPosition);
+
 private:
   struct InternalsType;
   std::shared_ptr<InternalsType> Internals;

vtkm/rendering/View.cxx

@@ -183,9 +183,6 @@ void View::RenderAnnotations()
 {
   if (this->RenderAnnotationsEnabled)
   {
-    this->SetupForScreenSpace();
-    this->RenderScreenAnnotations();
-
     this->GetCanvas().BeginTextRenderingBatch();
     for (auto& textAnnotation : this->Internal->TextAnnotations)
     {
@@ -203,6 +200,9 @@ void View::RenderAnnotations()
     {
       this->RenderWorldAnnotations();
     }
+
+    this->SetupForScreenSpace();
+    this->RenderScreenAnnotations();
   }
 }
 

vtkm/rendering/View.h

@@ -26,6 +26,9 @@ namespace vtkm
 {
 namespace rendering
 {
+const std::string RENDER_TIME_KEY = "RenderTime";
+const std::string COMPOSITE_TIME_KEY = "CompositeTime";
+const std::string TOTAL_TIME_KEY = "TotalTime";
 
 /// @brief The abstract class representing the view of a rendering scene.
 class VTKM_RENDERING_EXPORT View
@@ -133,6 +136,9 @@ public:
   VTKM_CONT
   void AddAdditionalAnnotation(std::function<void(void)> ann);
 
+  VTKM_CONT
+  std::unordered_map<std::string, vtkm::Float64> GetTimes() const { return this->Times; }
+
 protected:
   void SetupForWorldSpace(bool viewportClip = true);
 
@@ -143,6 +149,8 @@ protected:
   bool WorldAnnotationsEnabled = true;
   bool RenderAnnotationsEnabled = true;
 
+  std::unordered_map<std::string, vtkm::Float64> Times;
+
 private:
   std::unique_ptr<InternalData> Internal;
 };

vtkm/rendering/View3D.cxx

@@ -8,8 +8,16 @@
 //  PURPOSE.  See the above copyright notice for more information.
 //============================================================================
 
+#include <vtkm/cont/EnvironmentTracker.h>
+#include <vtkm/cont/Timer.h>
 #include <vtkm/rendering/View3D.h>
 
+#ifdef VTKM_ENABLE_MPI
+#include <mpi.h>
+#include <vtkm/thirdparty/diy/diy.h>
+#include <vtkm/thirdparty/diy/mpi-cast.h>
+#endif
+
 namespace vtkm
 {
 namespace rendering
@@ -36,20 +44,84 @@ View3D::View3D(const vtkm::rendering::Scene& scene,
 
 void View3D::Paint()
 {
+  this->Times.clear();
+
+  vtkm::cont::Timer totalTimer;
+  vtkm::cont::Timer renderTimer;
+  totalTimer.Start();
+  renderTimer.Start();
   this->GetCanvas().Clear();
-  this->RenderAnnotations();
   this->GetScene().Render(this->GetMapper(), this->GetCanvas(), this->GetCamera());
+  renderTimer.Stop();
+
+  vtkm::cont::Timer compositeTimer;
+  compositeTimer.Start();
+#ifdef VTKM_ENABLE_MPI
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+  if (comm.size() > 1)
+  {
+    this->Compositor.SetCompositeMode(vtkm::rendering::compositing::Compositor::Z_BUFFER_SURFACE);
+    this->Compositor.AddImage(this->GetCanvas());
+    auto result = this->Compositor.Composite();
+    //Rank 0 has the composited result, so put it into the Canvas.
+    if (comm.rank() == 0)
+    {
+      this->GetCanvas().CopyFrom(vtkm::cont::make_ArrayHandle(result.Pixels, vtkm::CopyFlag::Off),
+                                 vtkm::cont::make_ArrayHandle(result.Depths, vtkm::CopyFlag::Off));
+    }
+  }
+  this->RenderAnnotations();
+#endif
+  compositeTimer.Stop();
+  totalTimer.Stop();
+
+  this->Times[RENDER_TIME_KEY] = renderTimer.GetElapsedTime();
+  this->Times[COMPOSITE_TIME_KEY] = compositeTimer.GetElapsedTime();
+  this->Times[TOTAL_TIME_KEY] = totalTimer.GetElapsedTime();
 }
 
 void View3D::RenderScreenAnnotations()
 {
-  if (this->GetScene().GetNumberOfActors() > 0)
+  vtkm::Range scalarRange;
+
+  int numActors = this->GetScene().GetNumberOfActors();
+  if (numActors > 0)
+    scalarRange = this->GetScene().GetActor(0).GetScalarRange();
+
+  int totNumActors = numActors;
+
+  /*
+#ifdef VTKM_ENABLE_MPI
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+
+  vtkm::Float64 minVal = scalarRange.Min, maxVal = scalarRange.Max;
+
+  MPI_Comm mpiComm = vtkmdiy::mpi::mpi_cast(comm.handle());
+  int totNumActors = 0;
+  vtkm::Float64 minVal_res = 0, maxVal_res = 0;
+  MPI_Reduce(&numActors, &totNumActors, 1, MPI_INT, MPI_SUM, 0, mpiComm);
+  MPI_Reduce(&minVal, &minVal_res, 1, MPI_DOUBLE, MPI_MIN, 0, mpiComm);
+  MPI_Reduce(&maxVal, &maxVal_res, 1, MPI_DOUBLE, MPI_MAX, 0, mpiComm);
+  if (comm.rank() != 0)
+    return;
+
+  scalarRange.Min = minVal_res;
+  scalarRange.Max = maxVal_res;
+#endif
+
+  std::cout<<"totNumActors= "<<totNumActors<<" range= "<<scalarRange<<std::endl;
+
+  //DRP
+  //This assumes that rank 0 has an actor!!
+  */
+
+  if (totNumActors > 0)
   {
     this->GetCanvas().BeginTextRenderingBatch();
     this->GetWorldAnnotator().BeginLineRenderingBatch();
     //this->ColorBarAnnotation.SetAxisColor(vtkm::rendering::Color(1,1,1));
     this->ColorBarAnnotation.SetFieldName(this->GetScene().GetActor(0).GetScalarField().GetName());
-    this->ColorBarAnnotation.SetRange(this->GetScene().GetActor(0).GetScalarRange(), 5);
+    this->ColorBarAnnotation.SetRange(scalarRange, 5);
     this->ColorBarAnnotation.SetColorTable(this->GetScene().GetActor(0).GetColorTable());
     this->ColorBarAnnotation.Render(
       this->GetCamera(), this->GetWorldAnnotator(), this->GetCanvas());
@@ -60,8 +132,39 @@ void View3D::RenderScreenAnnotations()
 
 void View3D::RenderWorldAnnotations()
 {
-  this->GetCanvas().BeginTextRenderingBatch();
   vtkm::Bounds bounds = this->GetScene().GetSpatialBounds();
+
+#ifdef VTKM_ENABLE_MPI
+  //For parallel, get the collective bounds.
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+
+  vtkm::Float64 mins[3], maxs[3], mins_res[3], maxs_res[3];
+  mins[0] = bounds.X.Min;
+  mins[1] = bounds.Y.Min;
+  mins[2] = bounds.Z.Min;
+  maxs[0] = bounds.X.Max;
+  maxs[1] = bounds.Y.Max;
+  maxs[2] = bounds.Z.Max;
+
+  //DRP
+  //what if a scene has NO actors??
+
+  MPI_Comm mpiComm = vtkmdiy::mpi::mpi_cast(comm.handle());
+  MPI_Reduce(mins, mins_res, 3, MPI_DOUBLE, MPI_MIN, 0, mpiComm);
+  MPI_Reduce(maxs, maxs_res, 3, MPI_DOUBLE, MPI_MAX, 0, mpiComm);
+  if (comm.rank() != 0)
+    return;
+
+  bounds.X.Min = mins_res[0];
+  bounds.Y.Min = mins_res[1];
+  bounds.Z.Min = mins_res[2];
+  bounds.X.Max = maxs_res[0];
+  bounds.Y.Max = maxs_res[1];
+  bounds.Z.Max = maxs_res[2];
+
+#endif
+
+  this->GetCanvas().BeginTextRenderingBatch();
   vtkm::Float64 xmin = bounds.X.Min, xmax = bounds.X.Max;
   vtkm::Float64 ymin = bounds.Y.Min, ymax = bounds.Y.Max;
   vtkm::Float64 zmin = bounds.Z.Min, zmax = bounds.Z.Max;
@@ -70,7 +173,7 @@ void View3D::RenderWorldAnnotations()
 
   this->GetWorldAnnotator().BeginLineRenderingBatch();
   this->BoxAnnotation.SetColor(Color(.5f, .5f, .5f));
-  this->BoxAnnotation.SetExtents(this->GetScene().GetSpatialBounds());
+  this->BoxAnnotation.SetExtents(bounds);
   this->BoxAnnotation.Render(this->GetCamera(), this->GetWorldAnnotator());
   this->GetWorldAnnotator().EndLineRenderingBatch();
 

vtkm/rendering/View3D.h

@@ -16,6 +16,10 @@
 #include <vtkm/rendering/BoundingBoxAnnotation.h>
 #include <vtkm/rendering/ColorBarAnnotation.h>
 
+#ifdef VTKM_ENABLE_MPI
+#include <vtkm/rendering/compositing/Compositor.h>
+#endif
+
 namespace vtkm
 {
 namespace rendering
@@ -52,6 +56,10 @@ private:
   vtkm::rendering::AxisAnnotation3D YAxisAnnotation;
   vtkm::rendering::AxisAnnotation3D ZAxisAnnotation;
   vtkm::rendering::ColorBarAnnotation ColorBarAnnotation;
+
+#ifdef VTKM_ENABLE_MPI
+  vtkm::rendering::compositing::Compositor Compositor;
+#endif
 };
 }
 } // namespace vtkm::rendering

vtkm/rendering/compositing/AbsorptionPartial.h

@@ -0,0 +1,80 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtkm_render_compositing_Absorbtion_Partial_h
+#define vtkm_render_compositing_Absorbtion_Partial_h
+
+#include <assert.h>
+#include <vector>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+template <typename FloatType>
+struct AbsorptionPartial
+{
+  typedef FloatType ValueType;
+  int m_pixel_id;
+  double m_depth;
+  std::vector<FloatType> m_bins;
+
+  AbsorptionPartial()
+    : m_pixel_id(0)
+    , m_depth(0.f)
+  {
+  }
+
+  void print() {}
+
+  bool operator<(const AbsorptionPartial<FloatType>& other) const
+  {
+    //
+    // In absorption only we can blend the same
+    // pixel ids in any order
+    //
+    return m_pixel_id < other.m_pixel_id;
+  }
+
+  inline void blend(const AbsorptionPartial<FloatType>& other)
+  {
+    const int num_bins = static_cast<int>(m_bins.size());
+    assert(num_bins == (int)other.m_bins.size());
+    for (int i = 0; i < num_bins; ++i)
+    {
+      m_bins[i] *= other.m_bins[i];
+    }
+  }
+
+  static void composite_background(std::vector<AbsorptionPartial>& partials,
+                                   const std::vector<FloatType>& background)
+  {
+    const int size = static_cast<int>(partials.size());
+    AbsorptionPartial<FloatType> bg;
+    bg.m_bins = background;
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int i = 0; i < size; ++i)
+    {
+      partials[i].blend(bg);
+    }
+  }
+};
+
+}
+}
+} // vtkm::render::compositing
+
+
+#endif //vtkm_render_compositing_Absorbtion_Partial_h

vtkm/rendering/compositing/Compositor.cxx

@@ -0,0 +1,251 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#include <vtkm/cont/EnvironmentTracker.h>
+#include <vtkm/cont/ErrorBadValue.h>
+
+#include <algorithm>
+#include <vtkm/rendering/compositing/Compositor.h>
+#include <vtkm/rendering/compositing/ImageCompositor.h>
+
+#ifdef VTKM_ENABLE_MPI
+#include <mpi.h>
+#include <vtkm/rendering/compositing/DirectSendCompositor.h>
+#include <vtkm/rendering/compositing/RadixKCompositor.h>
+#include <vtkm/thirdparty/diy/mpi-cast.h>
+#endif
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+Compositor::Compositor()
+  : CompositingMode(Z_BUFFER_SURFACE)
+{
+}
+
+Compositor::~Compositor() {}
+
+void Compositor::SetCompositeMode(CompositeMode composite_mode)
+{
+  // assure we don't have mixed image types
+  assert(this->Images.size() == 0);
+  this->CompositingMode = composite_mode;
+}
+
+void Compositor::ClearImages()
+{
+  this->Images.clear();
+}
+
+void Compositor::AddImage(vtkm::rendering::Canvas& canvas)
+{
+  auto colors = &(canvas.GetColorBuffer().ReadPortal().GetArray()[0][0]);
+  auto depths = canvas.GetDepthBuffer().ReadPortal().GetArray();
+  vtkm::Id width = canvas.GetWidth();
+  vtkm::Id height = canvas.GetHeight();
+
+  // assert(this->CompositingMode != VIS_ORDER_BLEND);
+  assert(depths != NULL);
+  Image image;
+  if (this->Images.size() == 0)
+  {
+    this->Images.push_back(image);
+    this->Images[0].Init(colors, depths, width, height);
+    //this->Images[0].Save("first.png");
+  }
+  else if (this->CompositingMode == Z_BUFFER_SURFACE)
+  {
+    //
+    // Do local composite and keep a single image
+    //
+    image.Init(colors, depths, width, height);
+    vtkm::rendering::compositing::ImageCompositor compositor;
+    compositor.ZBufferComposite(this->Images[0], image);
+  }
+  else
+  {
+    const size_t image_index = this->Images.size();
+    this->Images.push_back(image);
+    this->Images[image_index].Init(colors, depths, width, height);
+  }
+}
+
+
+/*
+void Compositor::AddImage(const vtkm::cont::ArrayHandle<vtkm::Vec4<T>>& colors,
+                          const vtkm::cont::ArrayHandle<T>& depths,
+                          vtkm::Id width,
+                          vtkm::Id height)
+{
+  auto c = colors.WritePortal().GetArray();
+  auto d = depths.WritePortal().GetArray();
+  this->AddImage(c, d, width, height);
+}
+
+void Compositor::AddImage(const unsigned char* color_buffer,
+                          const float* depth_buffer,
+                          const int width,
+                          const int height)
+{
+  assert(this->CompositingMode != VIS_ORDER_BLEND);
+  assert(depth_buffer != NULL);
+  Image image;
+  if (this->Images.size() == 0)
+  {
+    this->Images.push_back(image);
+    this->Images[0].Init(color_buffer, depth_buffer, width, height);
+    //this->Images[0].Save("first.png");
+  }
+  else if (this->CompositingMode == Z_BUFFER_SURFACE)
+  {
+    //
+    // Do local composite and keep a single image
+    //
+    image.Init(color_buffer, depth_buffer, width, height);
+    vtkm::rendering::compositing::ImageCompositor compositor;
+    compositor.ZBufferComposite(this->Images[0], image);
+  }
+  else
+  {
+    const size_t image_index = this->Images.size();
+    this->Images.push_back(image);
+    this->Images[image_index].Init(color_buffer, depth_buffer, width, height);
+  }
+}
+
+void Compositor::AddImage(const float* color_buffer,
+                          const float* depth_buffer,
+                          const int width,
+                          const int height)
+{
+  assert(this->CompositingMode != VIS_ORDER_BLEND);
+  assert(depth_buffer != NULL);
+  Image image;
+  if (this->Images.size() == 0)
+  {
+    this->Images.push_back(image);
+    this->Images[0].Init(color_buffer, depth_buffer, width, height);
+  }
+  else if (this->CompositingMode == Z_BUFFER_SURFACE)
+  {
+    //
+    // Do local composite and keep a single image
+    //
+    image.Init(color_buffer, depth_buffer, width, height);
+
+    vtkm::rendering::compositing::ImageCompositor compositor;
+    compositor.ZBufferComposite(this->Images[0], image);
+  }
+  else
+  {
+    const size_t image_index = this->Images.size();
+    this->Images.push_back(image);
+    this->Images[image_index].Init(color_buffer, depth_buffer, width, height);
+  }
+}
+
+void Compositor::AddImage(const unsigned char* color_buffer,
+                          const float* depth_buffer,
+                          const int width,
+                          const int height,
+                          const int vis_order)
+{
+  assert(this->CompositingMode == VIS_ORDER_BLEND);
+  Image image;
+  const size_t image_index = this->Images.size();
+  this->Images.push_back(image);
+  this->Images[image_index].Init(color_buffer, depth_buffer, width, height, vis_order);
+}
+
+
+void Compositor::AddImage(const float* color_buffer,
+                          const float* depth_buffer,
+                          const int width,
+                          const int height,
+                          const int vis_order)
+{
+  assert(this->CompositingMode == VIS_ORDER_BLEND);
+  Image image;
+  const size_t image_index = this->Images.size();
+  this->Images.push_back(image);
+
+  this->Images[image_index].Init(color_buffer, depth_buffer, width, height, vis_order);
+}
+*/
+
+Image Compositor::Composite()
+{
+  assert(this->Images.size() != 0);
+
+  if (this->CompositingMode == Z_BUFFER_SURFACE)
+  {
+    CompositeZBufferSurface();
+  }
+  else if (this->CompositingMode == Z_BUFFER_BLEND)
+  {
+    CompositeZBufferBlend();
+  }
+  else if (this->CompositingMode == VIS_ORDER_BLEND)
+  {
+    CompositeVisOrder();
+  }
+  // Make this a param to avoid the copy?
+  return this->Images[0];
+}
+
+void Compositor::Cleanup() {}
+
+std::string Compositor::GetLogString()
+{
+  std::string res = m_log_stream.str();
+  m_log_stream.str("");
+  return res;
+}
+
+void Compositor::CompositeZBufferSurface()
+{
+  // nothing to do here in serial. Images were composited as
+  // they were added to the compositor
+#ifdef VTKM_ENABLE_MPI
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+
+  assert(this->Images.size() == 1);
+  RadixKCompositor compositor;
+  compositor.CompositeSurface(comm, this->Images[0]);
+  m_log_stream << compositor.GetTimingString();
+#endif
+}
+
+void Compositor::CompositeZBufferBlend()
+{
+  throw vtkm::cont::ErrorBadValue("Not implemented");
+}
+
+void Compositor::CompositeVisOrder()
+{
+
+#ifdef VTKM_ENABLE_MPI
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+  assert(this->Images.size() != 0);
+  vtkm::rendering::compositing::DirectSendCompositor compositor;
+  compositor.CompositeVolume(comm, this->Images);
+#else
+  vtkm::rendering::compositing::ImageCompositor compositor;
+  compositor.OrderedComposite(this->Images);
+#endif
+}
+
+}
+}
+} //namespace vtkm::rendering::compositing

vtkm/rendering/compositing/Compositor.h

@@ -0,0 +1,109 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_rendering_compositing_Compositor_h
+#define vtk_m_rendering_compositing_Compositor_h
+
+#include <vtkm/rendering/vtkm_rendering_export.h>
+
+#include <vtkm/rendering/Canvas.h>
+#include <vtkm/rendering/compositing/Image.h>
+
+#ifdef VTKM_ENABLE_MPI
+#include <mpi.h>
+#endif
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+class VTKM_RENDERING_EXPORT Compositor
+{
+public:
+  enum CompositeMode
+  {
+    Z_BUFFER_SURFACE, // zbuffer composite no transparency
+    Z_BUFFER_BLEND,   // zbuffer composite with transparency
+    VIS_ORDER_BLEND   // blend images in a specific order
+  };
+  Compositor();
+
+  virtual ~Compositor();
+
+  void SetCompositeMode(CompositeMode composite_mode);
+
+  void ClearImages();
+
+  void AddImage(vtkm::rendering::Canvas& canvas);
+
+  /*
+  void AddImage(const unsigned char* color_buffer,
+                const float* depth_buffer,
+                const int width,
+                const int height);
+
+  void AddImage(const float* color_buffer,
+                const float* depth_buffer,
+                const int width,
+                const int height);
+
+  void AddImage(const unsigned char* color_buffer,
+                const float* depth_buffer,
+                const int width,
+                const int height,
+                const int vis_order);
+
+  void AddImage(const float* color_buffer,
+                const float* depth_buffer,
+                const int width,
+                const int height,
+                const int vis_order);
+*/
+
+  Image Composite();
+
+  virtual void Cleanup();
+
+  std::string GetLogString();
+
+  unsigned char* ConvertBuffer(const float* buffer, const int size)
+  {
+    unsigned char* ubytes = new unsigned char[size];
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int i = 0; i < size; ++i)
+    {
+      ubytes[i] = static_cast<unsigned char>(buffer[i] * 255.f);
+    }
+
+    return ubytes;
+  }
+
+protected:
+  virtual void CompositeZBufferSurface();
+  virtual void CompositeZBufferBlend();
+  virtual void CompositeVisOrder();
+
+  std::stringstream m_log_stream;
+  CompositeMode CompositingMode;
+  std::vector<vtkm::rendering::compositing::Image> Images;
+};
+
+}
+}
+} //namespace vtkm::rendering::compositing
+
+
+#endif //vtk_m_rendering_compositing_Compositor_h

vtkm/rendering/compositing/DirectSendCompositor.cxx

@@ -0,0 +1,203 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#include <vtkm/rendering/compositing/DirectSendCompositor.h>
+#include <vtkm/rendering/compositing/ImageCompositor.h>
+
+#include <vtkm/rendering/compositing/vtkm_diy_collect.h>
+#include <vtkm/rendering/compositing/vtkm_diy_image_block.h>
+#include <vtkm/rendering/compositing/vtkm_diy_utils.h>
+
+/*
+#include <vtkh/compositing/MPICollect.hpp>
+#include <vtkh/compositing/vtkh_diy_collect.hpp>
+#include <vtkh/compositing/vtkh_diy_utils.hpp>
+
+#include <diy/master.hpp>
+#include <diy/mpi.hpp>
+#include <diy/partners/swap.hpp>
+#include <diy/reduce-operations.hpp>
+#include <diy/reduce.hpp>
+*/
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+namespace internal
+{
+struct Redistribute
+{
+  typedef vtkmdiy::RegularDecomposer<vtkmdiy::DiscreteBounds> Decomposer;
+  const vtkmdiy::RegularDecomposer<vtkmdiy::DiscreteBounds>& m_decomposer;
+  Redistribute(const Decomposer& decomposer)
+    : m_decomposer(decomposer)
+  {
+  }
+
+  void operator()(void* v_block, const vtkmdiy::ReduceProxy& proxy) const
+  {
+    MultiImageBlock* block = static_cast<MultiImageBlock*>(v_block);
+    //
+    // first round we have no incoming. Take the image we have,
+    // chop it up into pieces, and send it to the domain resposible
+    // for that portion
+    //
+    const int world_size = m_decomposer.nblocks;
+    const int local_images = block->m_images.size();
+    if (proxy.in_link().size() == 0)
+    {
+      std::map<vtkmdiy::BlockID, std::vector<Image>> outgoing;
+
+      for (int i = 0; i < world_size; ++i)
+      {
+        vtkmdiy::DiscreteBounds sub_image_bounds(3);
+        m_decomposer.fill_bounds(sub_image_bounds, i);
+        vtkm::Bounds vtkm_sub_bounds =
+          vtkm::rendering::compositing::DIYBoundsToVTKM(sub_image_bounds);
+
+        vtkmdiy::BlockID dest = proxy.out_link().target(i);
+        outgoing[dest].resize(local_images);
+
+        for (int img = 0; img < local_images; ++img)
+        {
+          outgoing[dest][img].SubsetFrom(block->m_images[img], vtkm_sub_bounds);
+        }
+      } //for
+
+      typename std::map<vtkmdiy::BlockID, std::vector<Image>>::iterator it;
+      for (it = outgoing.begin(); it != outgoing.end(); ++it)
+      {
+        proxy.enqueue(it->first, it->second);
+      }
+    } // if
+    else if (block->m_images.at(0).CompositeOrder != -1)
+    {
+      // blend images according to vis order
+      std::vector<Image> images;
+      for (int i = 0; i < proxy.in_link().size(); ++i)
+      {
+
+        std::vector<Image> incoming;
+        int gid = proxy.in_link().target(i).gid;
+        proxy.dequeue(gid, incoming);
+        const int in_size = incoming.size();
+        for (int img = 0; img < in_size; ++img)
+        {
+          images.emplace_back(incoming[img]);
+          //std::cout<<"rank "<<rank<<" rec "<<incoming[img].ToString()<<"\n";
+        }
+      } // for
+
+      ImageCompositor compositor;
+      compositor.OrderedComposite(images);
+
+      block->m_output.Swap(images[0]);
+    } // else if
+    else if (block->m_images.at(0).CompositeOrder == -1 &&
+             block->m_images.at(0).GetHasTransparency())
+    {
+      std::vector<Image> images;
+      for (int i = 0; i < proxy.in_link().size(); ++i)
+      {
+
+        std::vector<Image> incoming;
+        int gid = proxy.in_link().target(i).gid;
+        proxy.dequeue(gid, incoming);
+        const int in_size = incoming.size();
+        for (int img = 0; img < in_size; ++img)
+        {
+          images.emplace_back(incoming[img]);
+          //std::cout<<"rank "<<rank<<" rec "<<incoming[img].ToString()<<"\n";
+        }
+      } // for
+
+      //
+      // we have images with a depth buffer and transparency
+      //
+      ImageCompositor compositor;
+      compositor.ZBufferBlend(images);
+    }
+
+  } // operator
+};
+
+} //namespace internal
+
+DirectSendCompositor::DirectSendCompositor() {}
+
+DirectSendCompositor::~DirectSendCompositor() {}
+
+void DirectSendCompositor::CompositeVolume(vtkmdiy::mpi::communicator& diy_comm,
+                                           std::vector<Image>& images)
+{
+  vtkmdiy::DiscreteBounds global_bounds =
+    vtkm::rendering::compositing::VTKMBoundsToDIY(images.at(0).OrigBounds);
+
+  const int num_threads = 1;
+  const int num_blocks = diy_comm.size();
+  const int magic_k = 8;
+  Image sub_image;
+  //
+  // DIY does not seem to like being called with different block types
+  // so we isolate them within separate blocks
+  //
+  {
+    vtkmdiy::Master master(diy_comm, num_threads, -1, 0, [](void* b) {
+      ImageBlock<Image>* block = reinterpret_cast<ImageBlock<Image>*>(b);
+      delete block;
+    });
+
+    // create an assigner with one block per rank
+    vtkmdiy::ContiguousAssigner assigner(num_blocks, num_blocks);
+
+    AddMultiImageBlock create(master, images, sub_image);
+
+    const int dims = 2;
+    vtkmdiy::RegularDecomposer<vtkmdiy::DiscreteBounds> decomposer(dims, global_bounds, num_blocks);
+    decomposer.decompose(diy_comm.rank(), assigner, create);
+
+    vtkmdiy::all_to_all(master, assigner, internal::Redistribute(decomposer), magic_k);
+  }
+
+  {
+    vtkmdiy::Master master(diy_comm, num_threads, -1, 0, [](void* b) {
+      ImageBlock<Image>* block = reinterpret_cast<ImageBlock<Image>*>(b);
+      delete block;
+    });
+    vtkmdiy::ContiguousAssigner assigner(num_blocks, num_blocks);
+
+    const int dims = 2;
+    vtkmdiy::RegularDecomposer<vtkmdiy::DiscreteBounds> decomposer(dims, global_bounds, num_blocks);
+    AddImageBlock<Image> all_create(master, sub_image);
+    decomposer.decompose(diy_comm.rank(), assigner, all_create);
+    diy_comm.barrier();
+    //MPI_Barrier(diy_comm);
+
+    //MPICollect(sub_image,diy_comm);
+    vtkmdiy::all_to_all(master, assigner, CollectImages<Image>(decomposer), magic_k);
+  }
+
+  images.at(0).Swap(sub_image);
+}
+
+std::string DirectSendCompositor::GetTimingString()
+{
+  std::string res(m_timing_log.str());
+  m_timing_log.str("");
+  return res;
+}
+
+}
+}
+} //namespace vtkm::rendering::compositing

vtkm/rendering/compositing/DirectSendCompositor.h

@@ -0,0 +1,44 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_rendering_compositing_DirectSendCompositor_h
+#define vtk_m_rendering_compositing_DirectSendCompositor_h
+
+#include <vtkm/rendering/vtkm_rendering_export.h>
+
+#include <vtkm/rendering/compositing/Image.h>
+
+#include <vtkm/thirdparty/diy/diy.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+class VTKM_RENDERING_EXPORT DirectSendCompositor
+{
+public:
+  DirectSendCompositor();
+  ~DirectSendCompositor();
+  void CompositeVolume(vtkmdiy::mpi::communicator& diy_comm,
+                       std::vector<vtkm::rendering::compositing::Image>& images);
+  std::string GetTimingString();
+
+private:
+  std::stringstream m_timing_log;
+};
+
+}
+}
+} //namespace vtkm::rendering::compositing
+
+#endif //vtk_m_rendering_compositing_DirectSendCompositor_h

vtkm/rendering/compositing/EmissionPartial.h

@@ -0,0 +1,117 @@
+//============================================================================
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtkm_rendering_comositing_emmsion_partial_h
+#define vtkm_rendering_comositing_emmsion_partial_h
+
+#include <assert.h>
+#include <iostream>
+#include <vector>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+template <typename FloatType>
+struct EmissionPartial
+{
+  typedef FloatType ValueType;
+
+  int m_pixel_id;
+  double m_depth;
+  std::vector<FloatType> m_bins;
+  std::vector<FloatType> m_emission_bins;
+
+  EmissionPartial()
+    : m_pixel_id(0)
+    , m_depth(0.f)
+  {
+  }
+
+  void alter_bin(int bin, FloatType value)
+  {
+    m_bins[bin] = value;
+    m_emission_bins[bin] = value;
+  }
+
+  void print()
+  {
+    std::cout << "Partial id " << m_pixel_id << "\n";
+    std::cout << "Absorption : ";
+    for (int i = 0; i < m_bins.size(); ++i)
+    {
+      std::cout << m_bins[i] << " ";
+    }
+    std::cout << "\n";
+    std::cout << "Emission: ";
+    for (int i = 0; i < m_bins.size(); ++i)
+    {
+      std::cout << m_emission_bins[i] << " ";
+    }
+    std::cout << "\n";
+  }
+
+  bool operator<(const EmissionPartial<FloatType>& other) const
+  {
+    if (m_pixel_id != other.m_pixel_id)
+    {
+      return m_pixel_id < other.m_pixel_id;
+    }
+    else
+    {
+      return m_depth < other.m_depth;
+    }
+  }
+
+  inline void blend_absorption(const EmissionPartial<FloatType>& other)
+  {
+    const int num_bins = static_cast<int>(m_bins.size());
+    assert(num_bins == (int)other.m_bins.size());
+    for (int i = 0; i < num_bins; ++i)
+    {
+      m_bins[i] *= other.m_bins[i];
+    }
+  }
+
+  inline void blend_emission(EmissionPartial<FloatType>& other)
+  {
+    const int num_bins = static_cast<int>(m_bins.size());
+    assert(num_bins == (int)other.m_bins.size());
+    for (int i = 0; i < num_bins; ++i)
+    {
+      m_emission_bins[i] *= other.m_bins[i];
+    }
+  }
+
+  inline void add_emission(EmissionPartial<FloatType>& other)
+  {
+    const int num_bins = static_cast<int>(m_bins.size());
+    assert(num_bins == (int)other.m_bins.size());
+    for (int i = 0; i < num_bins; ++i)
+    {
+      m_emission_bins[i] += other.m_emission_bins[i];
+    }
+  }
+
+  static void composite_background(std::vector<EmissionPartial>& partials,
+                                   const std::vector<FloatType>& background)
+  {
+    //for(
+  }
+};
+
+}
+}
+} //vtkm::rendering::compositing
+
+
+#endif //vtkm_rendering_comositing_emmsion_partial_h

vtkm/rendering/compositing/Image.cxx

@@ -0,0 +1,35 @@
+// See License.txt
+
+#include <vtkm/rendering/compositing/Image.h>
+#include <vtkm/rendering/compositing/PNGEncoder.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+void Image::Save(const std::string& name, const std::vector<std::string>& comments)
+{
+  vtkm::rendering::compositing::PNGEncoder encoder;
+  encoder.Encode(&this->Pixels[0],
+                 this->Bounds.X.Max - this->Bounds.X.Min + 1,
+                 this->Bounds.Y.Max - this->Bounds.Y.Min + 1,
+                 comments);
+  encoder.Save(name);
+}
+
+void Image::Save(const std::string& name, const std::vector<std::string>& comments) const
+{
+  vtkm::rendering::compositing::PNGEncoder encoder;
+  encoder.Encode(&this->Pixels[0],
+                 this->Bounds.X.Max - this->Bounds.X.Min + 1,
+                 this->Bounds.Y.Max - this->Bounds.Y.Min + 1,
+                 comments);
+  encoder.Save(name);
+}
+
+}
+}
+} //namespace vtkm::rendering::compositing

vtkm/rendering/compositing/Image.h

@@ -0,0 +1,333 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_rendering_compositing_Image_h
+#define vtk_m_rendering_compositing_Image_h
+
+#include <vtkm/rendering/vtkm_rendering_export.h>
+
+#include <sstream>
+#include <vector>
+#include <vtkm/Bounds.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+struct VTKM_RENDERING_EXPORT Image
+{
+  // The image bounds are indicated by a grid starting at
+  // 1-width and 1-height. Actual width would be calculated
+  // Bounds.X.Max - Bounds.X.Min + 1
+  // 1024 - 1 + 1 = 1024
+  vtkm::Bounds OrigBounds;
+  vtkm::Bounds Bounds;
+  std::vector<unsigned char> Pixels;
+  std::vector<float> Depths;
+  int OrigRank;
+  bool HasTransparency;
+  int CompositeOrder;
+
+  Image()
+    : OrigRank(-1)
+    , HasTransparency(false)
+    , CompositeOrder(-1)
+  {
+  }
+
+
+  Image(const vtkm::Bounds& bounds)
+    : OrigBounds(bounds)
+    , Bounds(bounds)
+    , OrigRank(-1)
+    , HasTransparency(false)
+    , CompositeOrder(-1)
+
+  {
+    const int dx = bounds.X.Max - bounds.X.Min + 1;
+    const int dy = bounds.Y.Max - bounds.Y.Min + 1;
+    this->Pixels.resize(dx * dy * 4);
+    this->Depths.resize(dx * dy);
+  }
+
+  // init this image based on the original bounds
+  // of the other image
+  void InitOriginal(const Image& other)
+  {
+    this->OrigBounds = other.OrigBounds;
+    this->Bounds = other.OrigBounds;
+
+    const int dx = this->Bounds.X.Max - this->Bounds.X.Min + 1;
+    const int dy = this->Bounds.Y.Max - this->Bounds.Y.Min + 1;
+    this->Pixels.resize(dx * dy * 4);
+    this->Depths.resize(dx * dy);
+
+    this->OrigRank = -1;
+    this->HasTransparency = false;
+    this->CompositeOrder = -1;
+  }
+
+  int GetNumberOfPixels() const { return static_cast<int>(this->Pixels.size() / 4); }
+
+  void SetHasTransparency(bool has_transparency) { this->HasTransparency = has_transparency; }
+
+  bool GetHasTransparency() { return this->HasTransparency; }
+
+  void Init(const float* color_buffer,
+            const float* depth_buffer,
+            vtkm::Id width,
+            vtkm::Id height,
+            int composite_order = -1)
+  {
+    this->CompositeOrder = composite_order;
+    this->Bounds.X.Min = 1;
+    this->Bounds.Y.Min = 1;
+    this->Bounds.X.Max = width;
+    this->Bounds.Y.Max = height;
+    this->OrigBounds = this->Bounds;
+    const int size = width * height;
+    this->Pixels.resize(size * 4);
+    this->Depths.resize(size);
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int i = 0; i < size; ++i)
+    {
+      const int offset = i * 4;
+      this->Pixels[offset + 0] = static_cast<unsigned char>(color_buffer[offset + 0] * 255.f);
+      this->Pixels[offset + 1] = static_cast<unsigned char>(color_buffer[offset + 1] * 255.f);
+      this->Pixels[offset + 2] = static_cast<unsigned char>(color_buffer[offset + 2] * 255.f);
+      this->Pixels[offset + 3] = static_cast<unsigned char>(color_buffer[offset + 3] * 255.f);
+      float depth = depth_buffer[i];
+      //make sure we can do a single comparison on depth
+      //deal with negative depth values
+      //TODO: This may not be the best way
+      depth = depth < 0 ? abs(depth) : depth;
+      this->Depths[i] = depth;
+    }
+  }
+
+  void Init(const unsigned char* color_buffer,
+            const float* depth_buffer,
+            vtkm::Id width,
+            vtkm::Id height,
+            int composite_order = -1)
+  {
+    this->CompositeOrder = composite_order;
+    this->Bounds.X.Min = 1;
+    this->Bounds.Y.Min = 1;
+    this->Bounds.X.Max = width;
+    this->Bounds.Y.Max = height;
+    this->OrigBounds = this->Bounds;
+
+    const int size = width * height;
+    this->Pixels.resize(size * 4);
+    this->Depths.resize(size);
+
+    std::copy(color_buffer, color_buffer + size * 4, &this->Pixels[0]);
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int i = 0; i < size; ++i)
+    {
+      float depth = depth_buffer[i];
+      //make sure we can do a single comparison on depth
+      depth = depth < 0 ? 2.f : depth;
+      this->Depths[i] = depth;
+    } // for
+  }
+
+
+  void CompositeBackground(const float* color)
+  {
+
+    const int size = static_cast<int>(this->Pixels.size() / 4);
+    unsigned char bg_color[4];
+    for (int i = 0; i < 4; ++i)
+    {
+      bg_color[i] = static_cast<unsigned char>(color[i] * 255.f);
+    }
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int i = 0; i < size; ++i)
+    {
+      const int offset = i * 4;
+      unsigned int alpha = static_cast<unsigned int>(this->Pixels[offset + 3]);
+      const float opacity = (255 - alpha);
+      this->Pixels[offset + 0] += static_cast<unsigned char>(opacity * bg_color[0] / 255);
+      this->Pixels[offset + 1] += static_cast<unsigned char>(opacity * bg_color[1] / 255);
+      this->Pixels[offset + 2] += static_cast<unsigned char>(opacity * bg_color[2] / 255);
+      this->Pixels[offset + 3] += static_cast<unsigned char>(opacity * bg_color[3] / 255);
+    }
+  }
+  //
+  // Fill this image with a sub-region of another image
+  //
+  void SubsetFrom(const Image& image, const vtkm::Bounds& sub_region)
+  {
+    this->OrigBounds = image.OrigBounds;
+    this->Bounds = sub_region;
+    this->OrigRank = image.OrigRank;
+    this->CompositeOrder = image.CompositeOrder;
+
+    assert(sub_region.X.Min >= image.Bounds.X.Min);
+    assert(sub_region.Y.Min >= image.Bounds.Y.Min);
+    assert(sub_region.X.Max <= image.Bounds.X.Max);
+    assert(sub_region.Y.Max <= image.Bounds.Y.Max);
+
+    const int s_dx = this->Bounds.X.Max - this->Bounds.X.Min + 1;
+    const int s_dy = this->Bounds.Y.Max - this->Bounds.Y.Min + 1;
+
+    const int dx = image.Bounds.X.Max - image.Bounds.X.Min + 1;
+    //const int dy  = image.Bounds.Y.Max - image.Bounds.Y.Min + 1;
+
+    const int start_x = this->Bounds.X.Min - image.Bounds.X.Min;
+    const int start_y = this->Bounds.Y.Min - image.Bounds.Y.Min;
+    const int end_y = start_y + s_dy;
+
+    this->Pixels.resize(s_dx * s_dy * 4);
+    this->Depths.resize(s_dx * s_dy);
+
+
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int y = start_y; y < end_y; ++y)
+    {
+      const int copy_to = (y - start_y) * s_dx;
+      const int copy_from = y * dx + start_x;
+
+      std::copy(&image.Pixels[copy_from * 4],
+                &image.Pixels[copy_from * 4] + s_dx * 4,
+                &this->Pixels[copy_to * 4]);
+      std::copy(&image.Depths[copy_from], &image.Depths[copy_from] + s_dx, &this->Depths[copy_to]);
+    }
+  }
+
+  void Color(int color)
+  {
+    unsigned char c[4];
+    c[3] = 255;
+
+    c[0] = 0;
+    c[1] = 0;
+    c[2] = 0;
+    int index = color % 3;
+    c[index] = 255 - color * 11;
+    ;
+    const int size = static_cast<int>(this->Pixels.size());
+    for (int i = 0; i < size; ++i)
+    {
+      float d = this->Depths[i / 4];
+      if (d > 0 && d < 1)
+      {
+        this->Pixels[i] = c[i % 4];
+      }
+      else
+      {
+        this->Pixels[i] = 155;
+      }
+    }
+  }
+  //
+  // Fills the passed in image with the contents of this image
+  //
+  void SubsetTo(Image& image) const
+  {
+    image.CompositeOrder = this->CompositeOrder;
+    assert(this->Bounds.X.Min >= image.Bounds.X.Min);
+    assert(this->Bounds.Y.Min >= image.Bounds.Y.Min);
+    assert(this->Bounds.X.Max <= image.Bounds.X.Max);
+    assert(this->Bounds.Y.Max <= image.Bounds.Y.Max);
+
+    const int s_dx = this->Bounds.X.Max - this->Bounds.X.Min + 1;
+    const int s_dy = this->Bounds.Y.Max - this->Bounds.Y.Min + 1;
+
+    const int dx = image.Bounds.X.Max - image.Bounds.X.Min + 1;
+    //const int dy  = image.Bounds.Y.Max - image.Bounds.Y.Min + 1;
+
+    const int start_x = this->Bounds.X.Min - image.Bounds.X.Min;
+    const int start_y = this->Bounds.Y.Min - image.Bounds.Y.Min;
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int y = 0; y < s_dy; ++y)
+    {
+      const int copy_to = (y + start_y) * dx + start_x;
+      const int copy_from = y * s_dx;
+
+      std::copy(&this->Pixels[copy_from * 4],
+                &this->Pixels[copy_from * 4] + s_dx * 4,
+                &image.Pixels[copy_to * 4]);
+
+      std::copy(&this->Depths[copy_from], &this->Depths[copy_from] + s_dx, &image.Depths[copy_to]);
+    }
+  }
+
+  void Swap(Image& other)
+  {
+    vtkm::Bounds orig = this->OrigBounds;
+    vtkm::Bounds bounds = this->Bounds;
+
+    this->OrigBounds = other.OrigBounds;
+    this->Bounds = other.Bounds;
+
+    other.OrigBounds = orig;
+    other.Bounds = bounds;
+
+    this->Pixels.swap(other.Pixels);
+    this->Depths.swap(other.Depths);
+  }
+
+  void Clear()
+  {
+    vtkm::Bounds empty;
+    this->OrigBounds = empty;
+    this->Bounds = empty;
+    this->Pixels.clear();
+    this->Depths.clear();
+  }
+
+  std::string ToString() const
+  {
+    std::stringstream ss;
+    ss << "Total size pixels " << (int)this->Pixels.size() / 4;
+    ss << " tile dims: {" << this->Bounds.X.Min << "," << this->Bounds.Y.Min << "} - ";
+    ss << "{" << this->Bounds.X.Max << "," << this->Bounds.Y.Max << "}\n";
+    ;
+    return ss.str();
+  }
+
+  void Save(const std::string& name, const std::vector<std::string>& comments) const;
+  void Save(const std::string& name, const std::vector<std::string>& comments);
+};
+
+struct CompositeOrderSort
+{
+  inline bool operator()(const Image& lhs, const Image& rhs) const
+  {
+    return lhs.CompositeOrder < rhs.CompositeOrder;
+  }
+};
+
+}
+}
+} //namespace vtkm::rendering::compositing
+
+#endif //vtk_m_rendering_compositing_Image_h

vtkm/rendering/compositing/ImageCompositor.h

@@ -0,0 +1,217 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_rendering_compositing_ImageCompositor_h
+#define vtk_m_rendering_compositing_ImageCompositor_h
+
+#include <vtkm/rendering/vtkm_rendering_export.h>
+
+#include <algorithm>
+#include <vtkm/rendering/compositing/Image.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+class VTKM_RENDERING_EXPORT ImageCompositor
+{
+public:
+  void Blend(vtkm::rendering::compositing::Image& front, vtkm::rendering::compositing::Image& back)
+  {
+    assert(front.Bounds.X.Min == back.Bounds.X.Min);
+    assert(front.Bounds.Y.Min == back.Bounds.Y.Min);
+    assert(front.Bounds.X.Max == back.Bounds.X.Max);
+    assert(front.Bounds.Y.Max == back.Bounds.Y.Max);
+    const int size = static_cast<int>(front.Pixels.size() / 4);
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int i = 0; i < size; ++i)
+    {
+      const int offset = i * 4;
+      unsigned int alpha = front.Pixels[offset + 3];
+      const unsigned int opacity = 255 - alpha;
+
+      front.Pixels[offset + 0] +=
+        static_cast<unsigned char>(opacity * back.Pixels[offset + 0] / 255);
+      front.Pixels[offset + 1] +=
+        static_cast<unsigned char>(opacity * back.Pixels[offset + 1] / 255);
+      front.Pixels[offset + 2] +=
+        static_cast<unsigned char>(opacity * back.Pixels[offset + 2] / 255);
+      front.Pixels[offset + 3] +=
+        static_cast<unsigned char>(opacity * back.Pixels[offset + 3] / 255);
+
+      float d1 = std::min(front.Depths[i], 1.001f);
+      float d2 = std::min(back.Depths[i], 1.001f);
+      float depth = std::min(d1, d2);
+      front.Depths[i] = depth;
+    }
+  }
+
+  void ZBufferComposite(vtkm::rendering::compositing::Image& front,
+                        const vtkm::rendering::compositing::Image& image)
+  {
+    assert(front.Depths.size() == front.Pixels.size() / 4);
+    assert(front.Bounds.X.Min == image.Bounds.X.Min);
+    assert(front.Bounds.Y.Min == image.Bounds.Y.Min);
+    assert(front.Bounds.X.Max == image.Bounds.X.Max);
+    assert(front.Bounds.Y.Max == image.Bounds.Y.Max);
+
+    const int size = static_cast<int>(front.Depths.size());
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int i = 0; i < size; ++i)
+    {
+      const float depth = image.Depths[i];
+      if (depth > 1.f || front.Depths[i] < depth)
+      {
+        continue;
+      }
+      const int offset = i * 4;
+      front.Depths[i] = abs(depth);
+      front.Pixels[offset + 0] = image.Pixels[offset + 0];
+      front.Pixels[offset + 1] = image.Pixels[offset + 1];
+      front.Pixels[offset + 2] = image.Pixels[offset + 2];
+      front.Pixels[offset + 3] = image.Pixels[offset + 3];
+    }
+  }
+
+  void OrderedComposite(std::vector<vtkm::rendering::compositing::Image>& images)
+  {
+    const int total_images = images.size();
+    std::sort(images.begin(), images.end(), CompositeOrderSort());
+    for (int i = 1; i < total_images; ++i)
+    {
+      Blend(images[0], images[i]);
+    }
+  }
+
+  void ZBufferComposite(std::vector<vtkm::rendering::compositing::Image>& images)
+  {
+    const int total_images = images.size();
+    for (int i = 1; i < total_images; ++i)
+    {
+      ZBufferComposite(images[0], images[i]);
+    }
+  }
+
+  struct Pixel
+  {
+    unsigned char Color[4];
+    float Depth;
+    int PixelId; // local (sub-image) pixels id
+
+    bool operator<(const Pixel& other) const
+    {
+      if (this->PixelId != other.PixelId)
+      {
+        return this->PixelId < other.PixelId;
+      }
+      else
+      {
+        return this->Depth < other.Depth;
+      }
+    }
+  };
+
+  void CombineImages(const std::vector<vtkm::rendering::compositing::Image>& images,
+                     std::vector<Pixel>& pixels)
+  {
+
+    const int num_images = static_cast<int>(images.size());
+    for (int i = 0; i < num_images; ++i)
+    {
+      //
+      //  Extract the partial composites into a contiguous array
+      //
+
+      const int image_size = images[i].GetNumberOfPixels();
+      const int offset = i * image_size;
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+      for (int j = 0; j < image_size; ++j)
+      {
+        const int image_offset = j * 4;
+        pixels[offset + j].Color[0] = images[i].Pixels[image_offset + 0];
+        pixels[offset + j].Color[1] = images[i].Pixels[image_offset + 1];
+        pixels[offset + j].Color[2] = images[i].Pixels[image_offset + 2];
+        pixels[offset + j].Color[3] = images[i].Pixels[image_offset + 3];
+        pixels[offset + j].Depth = images[i].Depths[j];
+        pixels[offset + j].PixelId = j;
+      } // for pixels
+    }   // for images
+  }
+
+  void ZBufferBlend(std::vector<vtkm::rendering::compositing::Image>& images)
+  {
+    const int image_pixels = images[0].GetNumberOfPixels();
+    const int num_images = static_cast<int>(images.size());
+    std::vector<Pixel> pixels;
+    CombineImages(images, pixels);
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int i = 0; i < image_pixels; ++i)
+    {
+      const int begin = image_pixels * i;
+      const int end = image_pixels * i - 1;
+      std::sort(pixels.begin() + begin, pixels.begin() + end);
+    }
+
+    // check to see if that worked
+    int pixel_id_0 = pixels[0].PixelId;
+    for (int i = 1; i < num_images; ++i)
+    {
+      assert(pixel_id_0 == pixels[i].PixelId);
+    }
+
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int i = 0; i < image_pixels; ++i)
+    {
+      const int index = i * num_images;
+      Pixel pixel = pixels[index];
+      for (int j = 1; j < num_images; ++j)
+      {
+        if (pixel.Color[3] == 255 || pixel.Depth > 1.f)
+        {
+          break;
+        }
+        unsigned int alpha = pixel.Color[3];
+        const unsigned int opacity = 255 - alpha;
+        pixel.Color[0] += static_cast<unsigned char>(opacity * pixels[index + j].Color[0] / 255);
+        pixel.Color[1] += static_cast<unsigned char>(opacity * pixels[index + j].Color[1] / 255);
+        pixel.Color[2] += static_cast<unsigned char>(opacity * pixels[index + j].Color[2] / 255);
+        pixel.Color[3] += static_cast<unsigned char>(opacity * pixels[index + j].Color[3] / 255);
+        pixel.Depth = pixels[index + j].Depth;
+      } // for each image
+      images[0].Pixels[i * 4 + 0] = pixel.Color[0];
+      images[0].Pixels[i * 4 + 1] = pixel.Color[1];
+      images[0].Pixels[i * 4 + 2] = pixel.Color[2];
+      images[0].Pixels[i * 4 + 3] = pixel.Color[3];
+      images[0].Depths[i] = pixel.Depth;
+    } // for each pixel
+  }
+};
+
+}
+}
+} //namespace vtkm::rendering::compositing
+
+#endif //vtk_m_rendering_compositing_ImageComposititing_h

vtkm/rendering/compositing/PNGEncoder.cxx

@@ -0,0 +1,246 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#include <vtkm/internal/Configure.h>
+#include <vtkm/rendering/compositing/PNGEncoder.h>
+
+#include <iostream>
+#include <stdlib.h>
+
+VTKM_THIRDPARTY_PRE_INCLUDE
+#include <vtkm/thirdparty/lodepng/vtkmlodepng/lodepng.h>
+VTKM_THIRDPARTY_POST_INCLUDE
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+PNGEncoder::PNGEncoder()
+  : m_buffer(NULL)
+  , m_buffer_size(0)
+{
+}
+
+PNGEncoder::~PNGEncoder()
+{
+  Cleanup();
+}
+
+void PNGEncoder::Encode(const unsigned char* rgba_in, const int width, const int height)
+{
+  Cleanup();
+
+  // upside down relative to what lodepng wants
+  unsigned char* rgba_flip = new unsigned char[width * height * 4];
+
+  for (int y = 0; y < height; ++y)
+  {
+    memcpy(&(rgba_flip[y * width * 4]), &(rgba_in[(height - y - 1) * width * 4]), width * 4);
+  }
+
+  vtkm::png::LodePNGState state;
+  vtkm::png::lodepng_state_init(&state);
+  // use less aggressive compression
+  state.encoder.zlibsettings.btype = 2;
+  state.encoder.zlibsettings.use_lz77 = 0;
+
+  unsigned error = lodepng_encode(&m_buffer, &m_buffer_size, &rgba_flip[0], width, height, &state);
+  delete[] rgba_flip;
+
+  if (error)
+  {
+    std::cerr << "lodepng_encode_memory failed\n";
+  }
+}
+
+void PNGEncoder::Encode(const float* rgba_in, const int width, const int height)
+{
+  Cleanup();
+
+  // upside down relative to what lodepng wants
+  unsigned char* rgba_flip = new unsigned char[width * height * 4];
+
+
+  for (int x = 0; x < width; ++x)
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int y = 0; y < height; ++y)
+    {
+      int inOffset = (y * width + x) * 4;
+      int outOffset = ((height - y - 1) * width + x) * 4;
+      rgba_flip[outOffset + 0] = (unsigned char)(rgba_in[inOffset + 0] * 255.f);
+      rgba_flip[outOffset + 1] = (unsigned char)(rgba_in[inOffset + 1] * 255.f);
+      rgba_flip[outOffset + 2] = (unsigned char)(rgba_in[inOffset + 2] * 255.f);
+      rgba_flip[outOffset + 3] = (unsigned char)(rgba_in[inOffset + 3] * 255.f);
+    }
+
+  vtkm::png::LodePNGState state;
+  vtkm::png::lodepng_state_init(&state);
+  // use less aggressive compression
+  state.encoder.zlibsettings.btype = 2;
+  state.encoder.zlibsettings.use_lz77 = 0;
+
+  unsigned error = lodepng_encode(&m_buffer, &m_buffer_size, &rgba_flip[0], width, height, &state);
+  delete[] rgba_flip;
+
+  if (error)
+  {
+    std::cerr << "lodepng_encode_memory failed\n";
+  }
+}
+
+void PNGEncoder::Encode(const unsigned char* rgba_in,
+                        const int width,
+                        const int height,
+                        const std::vector<std::string>& comments)
+{
+  Cleanup();
+
+  // upside down relative to what lodepng wants
+  unsigned char* rgba_flip = new unsigned char[width * height * 4];
+
+  for (int y = 0; y < height; ++y)
+  {
+    memcpy(&(rgba_flip[y * width * 4]), &(rgba_in[(height - y - 1) * width * 4]), width * 4);
+  }
+
+  vtkm::png::LodePNGState state;
+  vtkm::png::lodepng_state_init(&state);
+  // use less aggressive compression
+  state.encoder.zlibsettings.btype = 2;
+  state.encoder.zlibsettings.use_lz77 = 0;
+  if (comments.size() % 2 != 0)
+  {
+    std::cerr << "PNGEncoder::Encode comments missing value for the last key.\n";
+    std::cerr << "Ignoring the last key.\n";
+  }
+  if (comments.size() > 1)
+  {
+    vtkm::png::lodepng_info_init(&state.info_png);
+    // Comments are in pairs with a key and a value, using
+    // comments.size()-1 ensures that we don't use the last
+    // comment if the length of the vector isn't a multiple of 2.
+    for (int i = 0; i < comments.size() - 1; i += 2)
+      vtkm::png::lodepng_add_text(&state.info_png, comments[i].c_str(), comments[i + 1].c_str());
+  }
+
+  unsigned error =
+    vtkm::png::lodepng_encode(&m_buffer, &m_buffer_size, &rgba_flip[0], width, height, &state);
+  delete[] rgba_flip;
+
+  if (error)
+  {
+    std::cerr << "lodepng_encode_memory failed\n";
+  }
+}
+
+void PNGEncoder::Encode(const float* rgba_in,
+                        const int width,
+                        const int height,
+                        const std::vector<std::string>& comments)
+{
+  Cleanup();
+
+  // upside down relative to what lodepng wants
+  unsigned char* rgba_flip = new unsigned char[width * height * 4];
+
+
+  for (int x = 0; x < width; ++x)
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int y = 0; y < height; ++y)
+    {
+      int inOffset = (y * width + x) * 4;
+      int outOffset = ((height - y - 1) * width + x) * 4;
+      rgba_flip[outOffset + 0] = (unsigned char)(rgba_in[inOffset + 0] * 255.f);
+      rgba_flip[outOffset + 1] = (unsigned char)(rgba_in[inOffset + 1] * 255.f);
+      rgba_flip[outOffset + 2] = (unsigned char)(rgba_in[inOffset + 2] * 255.f);
+      rgba_flip[outOffset + 3] = (unsigned char)(rgba_in[inOffset + 3] * 255.f);
+    }
+
+  vtkm::png::LodePNGState state;
+  vtkm::png::lodepng_state_init(&state);
+  // use less aggressive compression
+  state.encoder.zlibsettings.btype = 2;
+  state.encoder.zlibsettings.use_lz77 = 0;
+  if (comments.size() % 2 != 0)
+  {
+    std::cerr << "PNGEncoder::Encode comments missing value for the last key.\n";
+    std::cerr << "Ignoring the last key.\n";
+  }
+  if (comments.size() > 1)
+  {
+    vtkm::png::lodepng_info_init(&state.info_png);
+    // Comments are in pairs with a key and a value, using
+    // comments.size()-1 ensures that we don't use the last
+    // comment if the length of the vector isn't a multiple of 2.
+    for (int i = 0; i < comments.size() - 1; i += 2)
+      vtkm::png::lodepng_add_text(&state.info_png, comments[i].c_str(), comments[i + 1].c_str());
+  }
+
+  unsigned error =
+    vtkm::png::lodepng_encode(&m_buffer, &m_buffer_size, &rgba_flip[0], width, height, &state);
+  delete[] rgba_flip;
+
+  if (error)
+  {
+    std::cerr << "lodepng_encode_memory failed\n";
+  }
+}
+
+void PNGEncoder::Save(const std::string& filename)
+{
+  if (m_buffer == NULL)
+  {
+    std::cerr << "Save must be called after encode()\n";
+    /// we have a problem ...!
+    return;
+  }
+
+  unsigned error = vtkm::png::lodepng_save_file(m_buffer, m_buffer_size, filename.c_str());
+  if (error)
+  {
+    std::cerr << "Error saving PNG buffer to file: " << filename << "\n";
+  }
+}
+
+void* PNGEncoder::PngBuffer()
+{
+  return (void*)m_buffer;
+}
+
+size_t PNGEncoder::PngBufferSize()
+{
+  return m_buffer_size;
+}
+
+void PNGEncoder::Cleanup()
+{
+  if (m_buffer != NULL)
+  {
+    //lodepng_free(m_buffer);
+    // ^-- Not found even if LODEPNG_COMPILE_ALLOCATORS is defined?
+    // simply use "free"
+    free(m_buffer);
+    m_buffer = NULL;
+    m_buffer_size = 0;
+  }
+}
+
+}
+}
+} //namespace vtkm::rendering::compositing

vtkm/rendering/compositing/PNGEncoder.h

@@ -0,0 +1,57 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_rendering_compositing_PNGEncoder_h
+#define vtk_m_rendering_compositing_PNGEncoder_h
+
+#include <string>
+#include <vector>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+class PNGEncoder
+{
+public:
+  PNGEncoder();
+  ~PNGEncoder();
+
+  void Encode(const unsigned char* rgba_in, const int width, const int height);
+  void Encode(const float* rgba_in, const int width, const int height);
+  void Encode(const unsigned char* rgba_in,
+              const int width,
+              const int height,
+              const std::vector<std::string>& comments);
+  void Encode(const float* rgba_in,
+              const int width,
+              const int height,
+              const std::vector<std::string>& comments);
+  void Save(const std::string& filename);
+
+  void* PngBuffer();
+  size_t PngBufferSize();
+
+  void Cleanup();
+
+private:
+  unsigned char* m_buffer;
+  size_t m_buffer_size;
+};
+
+
+}
+}
+} //namespace vtkm::rendering::compositing
+
+#endif //vtk_m_rendering_compositing_PNGEncoder_h

vtkm/rendering/compositing/PartialCompositor.cxx

@@ -0,0 +1,552 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#include <vtkm/rendering/compositing/PartialCompositor.h>
+
+#include <algorithm>
+#include <assert.h>
+#include <limits>
+#include <vector>
+
+#ifdef VTKM_ENABLE_MPI
+#include <mpi.h>
+#include <vtkm/rendering/compositing/vtkm_diy_partial_collect.h>
+#include <vtkm/rendering/compositing/vtkm_diy_partial_redistribute.h>
+#endif
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+namespace detail
+{
+template <template <typename> class PartialType, typename FloatType>
+void BlendPartials(const int& total_segments,
+                   const int& total_partial_comps,
+                   std::vector<int>& pixel_work_ids,
+                   std::vector<PartialType<FloatType>>& partials,
+                   std::vector<PartialType<FloatType>>& output_partials,
+                   const int output_offset)
+{
+  //
+  // Perform the compositing and output the result in the output
+  //
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+  for (int i = 0; i < total_segments; ++i)
+  {
+    int current_index = pixel_work_ids[i];
+    PartialType<FloatType> result = partials[current_index];
+    ++current_index;
+    PartialType<FloatType> next = partials[current_index];
+    // TODO: we could just count the amount of work and make this a for loop(vectorize??)
+    while (result.m_pixel_id == next.m_pixel_id)
+    {
+      result.blend(next);
+      if (current_index + 1 >= total_partial_comps)
+      {
+        // we could break early for volumes,
+        // but blending past 1.0 alpha is no op.
+        break;
+      }
+      ++current_index;
+      next = partials[current_index];
+    }
+    output_partials[output_offset + i] = result;
+  }
+
+  //placeholder
+  //PartialType<FloatType>::composite_background(output_partials, background_values);
+}
+template <typename T>
+void BlendEmission(const int& total_segments,
+                   const int& total_partial_comps,
+                   std::vector<int>& pixel_work_ids,
+                   std::vector<EmissionPartial<T>>& partials,
+                   std::vector<EmissionPartial<T>>& output_partials,
+                   const int output_offset)
+{
+  //
+  // Perform the compositing and output the result in the output
+  // This code computes the optical depth (total absorption)
+  // along each rays path.
+  //
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+  for (int i = 0; i < total_segments; ++i)
+  {
+    int current_index = pixel_work_ids[i];
+    EmissionPartial<T> result = partials[current_index];
+    ++current_index;
+    EmissionPartial<T> next = partials[current_index];
+    // TODO: we could just count the amount of work and make this a for loop(vectorize??)
+    while (result.m_pixel_id == next.m_pixel_id)
+    {
+      result.blend_absorption(next);
+      if (current_index == total_partial_comps - 1)
+      {
+        break;
+      }
+      ++current_index;
+      next = partials[current_index];
+    }
+    output_partials[output_offset + i] = result;
+  }
+
+  //placeholder
+  //EmissionPartial::composite_background(output_partials);
+  // TODO: now blend source signature with output
+
+  //
+  //  Emission bins contain the amout of energy that leaves each
+  //  ray segment. To compute the amount of energy that reaches
+  //  the detector, we must multiply the segments emissed energy
+  //  by the optical depth of the remaining path to the detector.
+  //  To calculate the optical depth of the remaining path, we
+  //  do perform a reverse scan of absorption for each pixel id
+  //
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+  for (int i = 0; i < total_segments; ++i)
+  {
+    const int segment_start = pixel_work_ids[i];
+    int current_index = segment_start;
+    //
+    //  move forward to the end of the segment
+    //
+    while (partials[current_index].m_pixel_id == partials[current_index + 1].m_pixel_id)
+    {
+      ++current_index;
+      if (current_index == total_partial_comps - 1)
+      {
+        break;
+      }
+    }
+    //
+    // set the intensity emerging out of the last segment
+    //
+    output_partials[output_offset + i].m_emission_bins = partials[current_index].m_emission_bins;
+
+    //
+    // now move backwards accumulating absorption for each segment
+    // and then blending the intensity emerging from the previous
+    // segment.
+    //
+    current_index--;
+    while (current_index != segment_start - 1)
+    {
+      partials[current_index].blend_absorption(partials[current_index + 1]);
+      // mult this segments emission by the absorption in front
+      partials[current_index].blend_emission(partials[current_index + 1]);
+      // add remaining emissed engery to the output
+      output_partials[output_offset + i].add_emission(partials[current_index]);
+
+      --current_index;
+    }
+  }
+}
+template <>
+void BlendPartials<EmissionPartial, float>(const int& total_segments,
+                                           const int& total_partial_comps,
+                                           std::vector<int>& pixel_work_ids,
+                                           std::vector<EmissionPartial<float>>& partials,
+                                           std::vector<EmissionPartial<float>>& output_partials,
+                                           const int output_offset)
+{
+
+  BlendEmission(
+    total_segments, total_partial_comps, pixel_work_ids, partials, output_partials, output_offset);
+}
+
+template <>
+void BlendPartials<EmissionPartial, double>(const int& total_segments,
+                                            const int& total_partial_comps,
+                                            std::vector<int>& pixel_work_ids,
+                                            std::vector<EmissionPartial<double>>& partials,
+                                            std::vector<EmissionPartial<double>>& output_partials,
+                                            const int output_offset)
+{
+
+  BlendEmission(
+    total_segments, total_partial_comps, pixel_work_ids, partials, output_partials, output_offset);
+}
+
+} // namespace detail
+
+//--------------------------------------------------------------------------------------------
+template <typename PartialType>
+PartialCompositor<PartialType>::PartialCompositor()
+{
+}
+
+//--------------------------------------------------------------------------------------------
+
+template <typename PartialType>
+PartialCompositor<PartialType>::~PartialCompositor()
+{
+}
+
+//--------------------------------------------------------------------------------------------
+
+template <typename PartialType>
+void PartialCompositor<PartialType>::merge(const std::vector<std::vector<PartialType>>& in_partials,
+                                           std::vector<PartialType>& partials,
+                                           int& global_min_pixel,
+                                           int& global_max_pixel)
+{
+
+  int total_partial_comps = 0;
+  const int num_partial_images = static_cast<int>(in_partials.size());
+  int* offsets = new int[num_partial_images];
+  int* pixel_mins = new int[num_partial_images];
+  int* pixel_maxs = new int[num_partial_images];
+
+  for (int i = 0; i < num_partial_images; ++i)
+  {
+    offsets[i] = total_partial_comps;
+    total_partial_comps += in_partials[i].size();
+  }
+
+  partials.resize(total_partial_comps);
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+  for (int i = 0; i < num_partial_images; ++i)
+  {
+    //
+    //  Extract the partial composites into a contiguous array
+    //
+    std::copy(in_partials[i].begin(), in_partials[i].end(), partials.begin() + offsets[i]);
+  } // for each partial image
+
+  //
+  // Calculate the range of pixel ids
+  //
+  int max_pixel = std::numeric_limits<int>::min();
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for reduction(max : max_pixel)
+#endif
+  for (int i = 0; i < total_partial_comps; ++i)
+  {
+    int val = partials[i].m_pixel_id;
+    if (val > max_pixel)
+    {
+      max_pixel = val;
+    }
+  }
+
+  int min_pixel = std::numeric_limits<int>::max();
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for reduction(min : min_pixel)
+#endif
+  for (int i = 0; i < total_partial_comps; ++i)
+  {
+    int val = partials[i].m_pixel_id;
+    if (val < min_pixel)
+    {
+      min_pixel = val;
+    }
+  }
+
+  //
+  // determine the global pixel mins and maxs
+  //
+  global_min_pixel = min_pixel;
+  global_max_pixel = max_pixel;
+
+#ifdef VTKM_ENABLE_MPI
+  MPI_Comm comm_handle = MPI_Comm_f2c(m_mpi_comm_id);
+  int rank_min = global_min_pixel;
+  int rank_max = global_max_pixel;
+  int mpi_min;
+  int mpi_max;
+  MPI_Allreduce(&rank_min, &mpi_min, 1, MPI_INT, MPI_MIN, comm_handle);
+  MPI_Allreduce(&rank_max, &mpi_max, 1, MPI_INT, MPI_MAX, comm_handle);
+  global_min_pixel = mpi_min;
+  global_max_pixel = mpi_max;
+#endif
+
+  delete[] offsets;
+  delete[] pixel_mins;
+  delete[] pixel_maxs;
+}
+
+//--------------------------------------------------------------------------------------------
+template <typename PartialType>
+void PartialCompositor<PartialType>::composite_partials(std::vector<PartialType>& partials,
+                                                        std::vector<PartialType>& output_partials)
+{
+  const int total_partial_comps = partials.size();
+  if (total_partial_comps == 0)
+  {
+    output_partials = partials;
+    return;
+  }
+  //
+  // Sort the composites
+  //
+  std::sort(partials.begin(), partials.end());
+  //
+  // Find the number of unique pixel_ids with work
+  //
+  std::vector<unsigned char> work_flags;
+  std::vector<unsigned char> unique_flags;
+  work_flags.resize(total_partial_comps);
+  unique_flags.resize(total_partial_comps);
+  //
+  // just check the first and last entries manualy to reduce the
+  // loop complexity
+  //
+  if (partials[0].m_pixel_id == partials[1].m_pixel_id)
+  {
+    work_flags[0] = 1;
+    unique_flags[0] = 0;
+  }
+  else
+  {
+    work_flags[0] = 0;
+    unique_flags[0] = 1;
+  }
+  if (partials[total_partial_comps - 1].m_pixel_id != partials[total_partial_comps - 2].m_pixel_id)
+  {
+    unique_flags[total_partial_comps - 1] = 1;
+  }
+  else
+  {
+    unique_flags[total_partial_comps - 1] = 0;
+  }
+  const int n_minus_one = total_partial_comps - 1;
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+  for (int i = 1; i < n_minus_one; ++i)
+  {
+    unsigned char work_flag = 0;
+    unsigned char unique_flag = 0;
+    bool is_begining = false;
+    if (partials[i].m_pixel_id != partials[i - 1].m_pixel_id)
+    {
+      is_begining = true;
+    }
+
+    bool has_compositing_work = false;
+    if (partials[i].m_pixel_id == partials[i + 1].m_pixel_id)
+    {
+      has_compositing_work = true;
+    }
+    if (is_begining && has_compositing_work)
+    {
+      work_flag = 1;
+    }
+    if (is_begining && !has_compositing_work)
+    {
+      unique_flag = 1;
+    }
+
+    work_flags[i] = work_flag;
+    unique_flags[i] = unique_flag;
+  }
+  // count the number of of unique pixels
+  int total_segments = 0;
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for shared(work_flags) reduction(+ : total_segments)
+#endif
+  for (int i = 0; i < total_partial_comps; ++i)
+  {
+    total_segments += work_flags[i];
+  }
+
+  int total_unique_pixels = 0;
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for shared(unique_flags) reduction(+ : total_unique_pixels)
+#endif
+  for (int i = 0; i < total_partial_comps; ++i)
+  {
+    total_unique_pixels += unique_flags[i];
+  }
+
+  if (total_segments == 0)
+  {
+    //nothing to do
+  }
+
+  //
+  // find the pixel indexes that have compositing work
+  //
+  std::vector<int> pixel_work_ids;
+  pixel_work_ids.resize(total_segments);
+  int current_index = 0;
+  for (int i = 0; i < total_partial_comps; ++i)
+  {
+    if (work_flags[i] == 1)
+    {
+      pixel_work_ids[current_index] = i;
+      ++current_index;
+    }
+  }
+
+  //
+  // find the pixel indexes that have NO compositing work
+  //
+  std::vector<int> unique_ids;
+  unique_ids.resize(total_unique_pixels);
+  current_index = 0;
+  for (int i = 0; i < total_partial_comps; ++i)
+  {
+    if (unique_flags[i] == 1)
+    {
+      unique_ids[current_index] = i;
+      ++current_index;
+    }
+  }
+
+
+  const int total_output_pixels = total_unique_pixels + total_segments;
+
+  output_partials.resize(total_output_pixels);
+
+  //
+  // Gather the unique pixels into the output
+  //
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+  for (int i = 0; i < total_unique_pixels; ++i)
+  {
+    PartialType result = partials[unique_ids[i]];
+    output_partials[i] = result;
+  }
+
+  //
+  // perform compositing if there are more than
+  // one segment per ray
+  //
+  detail::BlendPartials(total_segments,
+                        total_partial_comps,
+                        pixel_work_ids,
+                        partials,
+                        output_partials,
+                        total_unique_pixels);
+}
+
+//--------------------------------------------------------------------------------------------
+
+template <typename PartialType>
+void PartialCompositor<PartialType>::composite(
+  std::vector<std::vector<PartialType>>& partial_images,
+  std::vector<PartialType>& output_partials)
+{
+  int global_partial_images = partial_images.size();
+#ifdef VTKM_ENABLE_MPI
+  MPI_Comm comm_handle = MPI_Comm_f2c(m_mpi_comm_id);
+  int local_partials = global_partial_images;
+  MPI_Allreduce(&local_partials, &global_partial_images, 1, MPI_INT, MPI_SUM, comm_handle);
+#endif
+
+#ifdef VTKM_ENABLE_MPI
+  // we could have no data, but it could exist elsewhere
+#endif
+
+  std::vector<PartialType> partials;
+  int global_min_pixel;
+  int global_max_pixel;
+
+  merge(partial_images, partials, global_min_pixel, global_max_pixel);
+
+  if (global_min_pixel > global_max_pixel)
+  {
+    // just bail
+    return;
+  }
+
+
+#ifdef VTKM_ENABLE_MPI
+  //
+  // Exchange partials with other ranks
+  //
+  redistribute(partials, comm_handle, global_min_pixel, global_max_pixel);
+  MPI_Barrier(comm_handle);
+#endif
+
+  const int total_partial_comps = partials.size();
+
+  //
+  // TODO: check to see if we have less than one
+  //
+  //assert(total_partial_comps > 1);
+
+  composite_partials(partials, output_partials);
+
+#ifdef VTKM_ENABLE_MPI
+  //
+  // Collect all of the distibuted pixels
+  //
+  collect(output_partials, comm_handle);
+  MPI_Barrier(comm_handle);
+#endif
+}
+
+template <typename PartialType>
+void PartialCompositor<PartialType>::set_background(std::vector<vtkm::Float32>& background_values)
+{
+  const size_t size = background_values.size();
+  m_background_values.resize(size);
+  for (size_t i = 0; i < size; ++i)
+  {
+    m_background_values[i] = background_values[i];
+  }
+}
+
+template <typename PartialType>
+void PartialCompositor<PartialType>::set_background(std::vector<vtkm::Float64>& background_values)
+{
+  const size_t size = background_values.size();
+  m_background_values.resize(size);
+  for (size_t i = 0; i < size; ++i)
+  {
+    m_background_values[i] = background_values[i];
+  }
+}
+
+template <typename PartialType>
+void PartialCompositor<PartialType>::set_comm_handle(int mpi_comm_id)
+{
+  m_mpi_comm_id = mpi_comm_id;
+}
+
+//Explicit function instantiations
+template class VTKM_RENDERING_EXPORT vtkm::rendering::compositing::PartialCompositor<
+  vtkm::rendering::compositing::VolumePartial<vtkm::Float32>>;
+template class VTKM_RENDERING_EXPORT vtkm::rendering::compositing::PartialCompositor<
+  vtkm::rendering::compositing::VolumePartial<vtkm::Float64>>;
+
+template class VTKM_RENDERING_EXPORT vtkm::rendering::compositing::PartialCompositor<
+  vtkm::rendering::compositing::AbsorptionPartial<vtkm::Float32>>;
+template class VTKM_RENDERING_EXPORT vtkm::rendering::compositing::PartialCompositor<
+  vtkm::rendering::compositing::AbsorptionPartial<vtkm::Float64>>;
+
+template class VTKM_RENDERING_EXPORT vtkm::rendering::compositing::PartialCompositor<
+  vtkm::rendering::compositing::EmissionPartial<vtkm::Float32>>;
+template class VTKM_RENDERING_EXPORT vtkm::rendering::compositing::PartialCompositor<
+  vtkm::rendering::compositing::EmissionPartial<vtkm::Float64>>;
+
+
+}
+}
+} //vtkm::rendering::compositing

vtkm/rendering/compositing/PartialCompositor.h

@@ -0,0 +1,60 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtkm_rendering_compositing_partial_compositor_h
+#define vtkm_rendering_compositing_partial_compositor_h
+
+#include <vtkm/rendering/vtkm_rendering_export.h>
+
+#include <iostream>
+#include <vector>
+#include <vtkm/Types.h>
+#include <vtkm/rendering/compositing/AbsorptionPartial.h>
+#include <vtkm/rendering/compositing/EmissionPartial.h>
+#include <vtkm/rendering/compositing/VolumePartial.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+template <typename PartialType>
+class VTKM_RENDERING_EXPORT PartialCompositor
+{
+public:
+  PartialCompositor();
+  ~PartialCompositor();
+  void composite(std::vector<std::vector<PartialType>>& partial_images,
+                 std::vector<PartialType>& output_partials);
+  void set_background(std::vector<vtkm::Float32>& background_values);
+  void set_background(std::vector<vtkm::Float64>& background_values);
+  void set_comm_handle(int mpi_comm_id);
+
+protected:
+  void merge(const std::vector<std::vector<PartialType>>& in_partials,
+             std::vector<PartialType>& partials,
+             int& global_min_pixel,
+             int& global_max_pixel);
+
+  void composite_partials(std::vector<PartialType>& partials,
+                          std::vector<PartialType>& output_partials);
+
+  std::vector<typename PartialType::ValueType> m_background_values;
+  int m_mpi_comm_id;
+};
+
+
+}
+}
+} //vtkm::rendering::compositing
+
+#endif //vtkm_rendering_compositing_partial_compositor_h

vtkm/rendering/compositing/PayloadCompositor.cxx

@@ -0,0 +1,80 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#include <vtkm/cont/EnvironmentTracker.h>
+#include <vtkm/rendering/compositing/PayloadCompositor.h>
+#include <vtkm/rendering/compositing/PayloadImageCompositor.h>
+#include <vtkm/rendering/compositing/RadixKCompositor.h>
+
+#include <algorithm>
+#include <assert.h>
+
+#ifdef VTKM_ENABLE_MPI
+#include <mpi.h>
+//#include <vtkh/vtkh.hpp>
+#include <vtkm/rendering/compositing/RadixKCompositor.h>
+//#include <diy/mpi.hpp>
+#endif
+
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+PayloadCompositor::PayloadCompositor() {}
+
+void PayloadCompositor::ClearImages()
+{
+  m_images.clear();
+}
+
+void PayloadCompositor::AddImage(vtkm::rendering::compositing::PayloadImage& image)
+{
+  assert(image.GetNumberOfPixels() != 0);
+
+  if (m_images.size() == 0)
+  {
+    m_images.push_back(image);
+  }
+  else
+  {
+    //
+    // Do local composite and keep a single image
+    //
+    vtkm::rendering::compositing::PayloadImageCompositor compositor;
+    compositor.ZBufferComposite(m_images[0], image);
+  }
+}
+
+vtkm::rendering::compositing::PayloadImage PayloadCompositor::Composite()
+{
+  assert(m_images.size() != 0);
+  // nothing to do here in serial. Images were composited as
+  // they were added to the compositor
+#ifdef VTKM_ENABLE_MPI
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+  //  vtkmdiy::mpi::communicator diy_comm;
+  //  diy_comm = vtkmdiy::mpi::communicator(MPI_Comm_f2c(GetMPICommHandle()));
+
+  assert(m_images.size() == 1);
+  vtkm::rendering::compositing::RadixKCompositor compositor;
+  compositor.CompositeSurface(comm, this->m_images[0]);
+#endif
+  // Make this a param to avoid the copy?
+  return m_images[0];
+}
+
+
+}
+}
+} // namespace vtkm:rendering::compositing

vtkm/rendering/compositing/PayloadCompositor.h

@@ -0,0 +1,45 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_rendering_compositing_PayloadCompositor_h
+#define vtk_m_rendering_compositing_PayloadCompositor_h
+
+#include <vtkm/rendering/vtkm_rendering_export.h>
+
+#include <vtkm/rendering/compositing/PayloadImage.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+class VTKM_RENDERING_EXPORT PayloadCompositor
+{
+public:
+  PayloadCompositor();
+
+  void ClearImages();
+
+  void AddImage(vtkm::rendering::compositing::PayloadImage& image);
+
+  vtkm::rendering::compositing::PayloadImage Composite();
+
+protected:
+  std::vector<vtkm::rendering::compositing::PayloadImage> m_images;
+};
+
+}
+}
+} // namespace vtkm:rendering::compositing
+
+
+#endif //vtk_m_rendering_compositing_PayloadCompositor_h

vtkm/rendering/compositing/PayloadImage.cxx

@@ -0,0 +1,33 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#include <vtkm/rendering/compositing/PNGEncoder.h>
+#include <vtkm/rendering/compositing/PayloadImage.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+void PayloadImage::Save(const std::string& name, const std::vector<std::string>& comments)
+{
+  PNGEncoder encoder;
+  encoder.Encode(&this->Payloads[0],
+                 this->Bounds.X.Max - this->Bounds.X.Min + 1,
+                 this->Bounds.Y.Max - this->Bounds.Y.Min + 1,
+                 comments);
+  encoder.Save(name);
+}
+
+}
+}
+} //namespace vtkm::rendering::compositing

vtkm/rendering/compositing/PayloadImage.h

@@ -0,0 +1,211 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_rendering_compositing_PayloadImage_h
+#define vtk_m_rendering_compositing_PayloadImage_h
+
+#include <vtkm/rendering/vtkm_rendering_export.h>
+
+#include <sstream>
+#include <vector>
+#include <vtkm/Bounds.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+struct VTKM_RENDERING_EXPORT PayloadImage
+{
+  // The image bounds are indicated by a grid starting at
+  // 1-width and 1-height. Actual width would be calculated
+  // Bounds.X.Max - Bounds.X.Min + 1
+  // 1024 - 1 + 1 = 1024
+  vtkm::Bounds OrigBounds;
+  vtkm::Bounds Bounds;
+  std::vector<unsigned char> Payloads;
+  std::vector<float> Depths;
+  int OrigRank;
+  int PayloadBytes; // Size of the payload in bytes
+  float DefaultValue;
+
+  PayloadImage() {}
+
+  PayloadImage(const vtkm::Bounds& bounds, const int payload_bytes)
+    : OrigBounds(bounds)
+    , Bounds(bounds)
+    , OrigRank(-1)
+    , PayloadBytes(payload_bytes)
+  {
+    DefaultValue = vtkm::Nan32();
+    const int dx = bounds.X.Max - bounds.X.Min + 1;
+    const int dy = bounds.Y.Max - bounds.Y.Min + 1;
+    this->Payloads.resize(dx * dy * this->PayloadBytes);
+    this->Depths.resize(dx * dy);
+  }
+
+  void InitOriginal(const PayloadImage& other)
+  {
+    this->OrigBounds = other.OrigBounds;
+    this->Bounds = other.OrigBounds;
+    this->PayloadBytes = other.PayloadBytes;
+    this->DefaultValue = other.DefaultValue;
+
+    const int dx = this->Bounds.X.Max - this->Bounds.X.Min + 1;
+    const int dy = this->Bounds.Y.Max - this->Bounds.Y.Min + 1;
+    this->Payloads.resize(dx * dy * this->PayloadBytes);
+    this->Depths.resize(dx * dy);
+
+    this->OrigRank = -1;
+  }
+
+  int GetNumberOfPixels() const { return static_cast<int>(this->Depths.size()); }
+
+  void Init(const unsigned char* payload_buffer, const float* depth_buffer, int width, int height)
+  {
+    this->Bounds.X.Min = 1;
+    this->Bounds.Y.Min = 1;
+    this->Bounds.X.Max = width;
+    this->Bounds.Y.Max = height;
+    this->OrigBounds = this->Bounds;
+    const int size = width * height;
+    this->Payloads.resize(size * this->PayloadBytes);
+    this->Depths.resize(size);
+
+    std::copy(payload_buffer, payload_buffer + size * this->PayloadBytes, &this->Payloads[0]);
+
+    std::copy(depth_buffer, depth_buffer + size, &this->Depths[0]);
+  }
+
+  //
+  // Fill this image with a sub-region of another image
+  //
+  void SubsetFrom(const PayloadImage& image, const vtkm::Bounds& sub_region)
+  {
+    this->OrigBounds = image.OrigBounds;
+    this->Bounds = sub_region;
+    this->OrigRank = image.OrigRank;
+    this->PayloadBytes = image.PayloadBytes;
+
+    assert(sub_region.X.Min >= image.Bounds.X.Min);
+    assert(sub_region.Y.Min >= image.Bounds.Y.Min);
+    assert(sub_region.X.Max <= image.Bounds.X.Max);
+    assert(sub_region.Y.Max <= image.Bounds.Y.Max);
+
+    const int s_dx = this->Bounds.X.Max - this->Bounds.X.Min + 1;
+    const int s_dy = this->Bounds.Y.Max - this->Bounds.Y.Min + 1;
+
+    const int dx = image.Bounds.X.Max - image.Bounds.X.Min + 1;
+    //const int dy  = image.Bounds.Y.Max - image.Bounds.Y.Min + 1;
+
+    const int start_x = this->Bounds.X.Min - image.Bounds.X.Min;
+    const int start_y = this->Bounds.Y.Min - image.Bounds.Y.Min;
+    const int end_y = start_y + s_dy;
+
+    size_t buffer_size = s_dx * s_dy * this->PayloadBytes;
+
+    this->Payloads.resize(buffer_size);
+    this->Depths.resize(s_dx * s_dy);
+
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int y = start_y; y < end_y; ++y)
+    {
+      const int copy_to = (y - start_y) * s_dx;
+      const int copy_from = y * dx + start_x;
+
+      std::copy(&image.Payloads[copy_from * this->PayloadBytes],
+                &image.Payloads[copy_from * this->PayloadBytes] + s_dx * this->PayloadBytes,
+                &this->Payloads[copy_to * this->PayloadBytes]);
+      std::copy(&image.Depths[copy_from], &image.Depths[copy_from] + s_dx, &this->Depths[copy_to]);
+    }
+  }
+
+  //
+  // Fills the passed in image with the contents of this image
+  //
+  void SubsetTo(PayloadImage& image) const
+  {
+    assert(this->Bounds.X.Min >= image.Bounds.X.Min);
+    assert(this->Bounds.Y.Min >= image.Bounds.Y.Min);
+    assert(this->Bounds.X.Max <= image.Bounds.X.Max);
+    assert(this->Bounds.Y.Max <= image.Bounds.Y.Max);
+
+    const int s_dx = this->Bounds.X.Max - this->Bounds.X.Min + 1;
+    const int s_dy = this->Bounds.Y.Max - this->Bounds.Y.Min + 1;
+
+    const int dx = image.Bounds.X.Max - image.Bounds.X.Min + 1;
+    //const int dy  = image.Bounds.Y.Max - image.Bounds.Y.Min + 1;
+
+    const int start_x = this->Bounds.X.Min - image.Bounds.X.Min;
+    const int start_y = this->Bounds.Y.Min - image.Bounds.Y.Min;
+
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int y = 0; y < s_dy; ++y)
+    {
+      const int copy_to = (y + start_y) * dx + start_x;
+      const int copy_from = y * s_dx;
+
+      std::copy(&this->Payloads[copy_from * this->PayloadBytes],
+                &this->Payloads[copy_from * this->PayloadBytes] + s_dx * this->PayloadBytes,
+                &image.Payloads[copy_to * this->PayloadBytes]);
+
+      std::copy(&this->Depths[copy_from], &this->Depths[copy_from] + s_dx, &image.Depths[copy_to]);
+    }
+  }
+
+  void Swap(PayloadImage& other)
+  {
+    vtkm::Bounds orig = this->OrigBounds;
+    vtkm::Bounds bounds = this->Bounds;
+
+    this->OrigBounds = other.OrigBounds;
+    this->Bounds = other.Bounds;
+
+    other.OrigBounds = orig;
+    other.Bounds = bounds;
+
+    this->Payloads.swap(other.Payloads);
+    this->Depths.swap(other.Depths);
+  }
+
+  void Clear()
+  {
+    vtkm::Bounds empty;
+    this->OrigBounds = empty;
+    this->Bounds = empty;
+    this->Payloads.clear();
+    this->Depths.clear();
+  }
+
+  std::string ToString() const
+  {
+    std::stringstream ss;
+    ss << "Total size pixels " << (int)this->Depths.size();
+    ss << " tile dims: {" << this->Bounds.X.Min << "," << this->Bounds.Y.Min << "} - ";
+    ss << "{" << this->Bounds.X.Max << "," << this->Bounds.Y.Max << "}\n";
+    ;
+    return ss.str();
+  }
+
+  void Save(const std::string& name, const std::vector<std::string>& comments);
+};
+
+}
+}
+} //namespace vtkm::rendering::compositing
+
+#endif //vtk_m_rendering_compositing_PayloadImage_h

vtkm/rendering/compositing/PayloadImageCompositor.h

@@ -0,0 +1,74 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_rendering_compositing_PayloadImageCompositor_h
+#define vtk_m_rendering_compositing_PayloadImageCompositor_h
+
+#include <vtkm/rendering/vtkm_rendering_export.h>
+
+
+#include <algorithm>
+#include <cmath>
+#include <vtkm/rendering/compositing/PayloadImage.h>
+
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+class VTKM_RENDERING_EXPORT PayloadImageCompositor
+{
+public:
+  void ZBufferComposite(vtkm::rendering::compositing::PayloadImage& front,
+                        const vtkm::rendering::compositing::PayloadImage& image)
+  {
+    if (front.PayloadBytes != image.PayloadBytes)
+    {
+      std::cout << "very bad\n";
+    }
+    assert(front.Depths.size() == front.Payloads.size() / front.PayloadBytes);
+    assert(front.Bounds.X.Min == image.Bounds.X.Min);
+    assert(front.Bounds.Y.Min == image.Bounds.Y.Min);
+    assert(front.Bounds.X.Max == image.Bounds.X.Max);
+    assert(front.Bounds.Y.Max == image.Bounds.Y.Max);
+
+    const int size = static_cast<int>(front.Depths.size());
+    const bool nan_check = image.DefaultValue != image.DefaultValue;
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int i = 0; i < size; ++i)
+    {
+      const float depth = image.Depths[i];
+      const float fdepth = front.Depths[i];
+      // this should handle NaNs correctly
+      const bool take_back = fmin(depth, fdepth) == depth;
+
+      if (take_back)
+      {
+        const int offset = i * 4;
+        front.Depths[i] = depth;
+        const size_t p_offset = i * front.PayloadBytes;
+        std::copy(&image.Payloads[p_offset],
+                  &image.Payloads[p_offset] + front.PayloadBytes,
+                  &front.Payloads[p_offset]);
+      }
+    }
+  }
+};
+
+}
+}
+} //namespace vtkm::rendering::compositing
+
+#endif //vtk_m_rendering_compositing_PayloadImageCompositor_h

vtkm/rendering/compositing/RadixKCompositor.cxx

@@ -0,0 +1,210 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#include <vtkm/rendering/compositing/Image.h>
+#include <vtkm/rendering/compositing/ImageCompositor.h>
+#include <vtkm/rendering/compositing/PayloadImage.h>
+#include <vtkm/rendering/compositing/PayloadImageCompositor.h>
+//#include <vtkm/compositing/MPICollect.hpp>
+#include <vtkm/rendering/compositing/RadixKCompositor.h>
+#include <vtkm/rendering/compositing/vtkm_diy_collect.h>
+#include <vtkm/rendering/compositing/vtkm_diy_utils.h>
+
+#include <vtkm/thirdparty/diy/diy.h>
+#include <vtkm/thirdparty/diy/master.h>
+//#include <vtkm/thirdparty/diy/mpi.h>
+#include <vtkm/thirdparty/diy/reduce-operations.h>
+#include <vtkm/thirdparty/diy/reduce.h>
+#include <vtkm/thirdparty/diy/swap.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+template <typename ImageType>
+void DepthComposite(ImageType& front, ImageType& back);
+
+template <>
+void DepthComposite<PayloadImage>(PayloadImage& front, PayloadImage& back)
+{
+  vtkm::rendering::compositing::PayloadImageCompositor compositor;
+  compositor.ZBufferComposite(front, back);
+}
+
+template <>
+void DepthComposite<vtkm::rendering::compositing::Image>(vtkm::rendering::compositing::Image& front,
+                                                         vtkm::rendering::compositing::Image& back)
+{
+  vtkm::rendering::compositing::ImageCompositor compositor;
+  compositor.ZBufferComposite(front, back);
+}
+
+template <typename ImageType>
+void reduce_images(void* b,
+                   const vtkmdiy::ReduceProxy& proxy,
+                   const vtkmdiy::RegularSwapPartners& partners)
+{
+  ImageBlock<ImageType>* block = reinterpret_cast<ImageBlock<ImageType>*>(b);
+  unsigned int round = proxy.round();
+  ImageType& image = block->m_image;
+  // count the number of incoming pixels
+  if (proxy.in_link().size() > 0)
+  {
+    for (int i = 0; i < proxy.in_link().size(); ++i)
+    {
+      int gid = proxy.in_link().target(i).gid;
+      if (gid == proxy.gid())
+      {
+        //skip revieving from self since we sent nothing
+        continue;
+      }
+      ImageType incoming;
+      proxy.dequeue(gid, incoming);
+      DepthComposite(image, incoming);
+    } // for in links
+  }
+
+  if (proxy.out_link().size() == 0)
+  {
+    return;
+  }
+  // do compositing?? intermediate stage?
+  const int group_size = proxy.out_link().size();
+  const int current_dim = partners.dim(round);
+
+  //create balanced set of ranges for current dim
+  vtkmdiy::DiscreteBounds image_bounds =
+    vtkm::rendering::compositing::VTKMBoundsToDIY(image.Bounds);
+  int range_length = image_bounds.max[current_dim] - image_bounds.min[current_dim];
+  int base_step = range_length / group_size;
+  int rem = range_length % group_size;
+  std::vector<int> bucket_sizes(group_size, base_step);
+  for (int i = 0; i < rem; ++i)
+  {
+    bucket_sizes[i]++;
+  }
+
+  int count = 0;
+  for (int i = 0; i < group_size; ++i)
+  {
+    count += bucket_sizes[i];
+  }
+  assert(count == range_length);
+
+  std::vector<vtkmdiy::DiscreteBounds> subset_bounds(
+    group_size, vtkm::rendering::compositing::VTKMBoundsToDIY(image.Bounds));
+  int min_pixel = image_bounds.min[current_dim];
+  for (int i = 0; i < group_size; ++i)
+  {
+    subset_bounds[i].min[current_dim] = min_pixel;
+    subset_bounds[i].max[current_dim] = min_pixel + bucket_sizes[i];
+    min_pixel += bucket_sizes[i];
+  }
+
+  //debug
+  if (group_size > 1)
+  {
+    for (int i = 1; i < group_size; ++i)
+    {
+      assert(subset_bounds[i - 1].max[current_dim] == subset_bounds[i].min[current_dim]);
+    }
+
+    assert(subset_bounds[0].min[current_dim] == image_bounds.min[current_dim]);
+    assert(subset_bounds[group_size - 1].max[current_dim] == image_bounds.max[current_dim]);
+  }
+
+  std::vector<ImageType> out_images(group_size);
+  for (int i = 0; i < group_size; ++i)
+  {
+    out_images[i].SubsetFrom(image,
+                             vtkm::rendering::compositing::DIYBoundsToVTKM(subset_bounds[i]));
+  } //for
+
+  for (int i = 0; i < group_size; ++i)
+  {
+    if (proxy.out_link().target(i).gid == proxy.gid())
+    {
+      image.Swap(out_images[i]);
+    }
+    else
+    {
+      proxy.enqueue(proxy.out_link().target(i), out_images[i]);
+    }
+  } //for
+
+} // reduce images
+
+RadixKCompositor::RadixKCompositor() {}
+
+RadixKCompositor::~RadixKCompositor() {}
+
+template <typename ImageType>
+void RadixKCompositor::CompositeImpl(vtkmdiy::mpi::communicator& diy_comm, ImageType& image)
+{
+  vtkmdiy::DiscreteBounds global_bounds =
+    vtkm::rendering::compositing::VTKMBoundsToDIY(image.OrigBounds);
+
+  // tells diy to use one thread
+  const int num_threads = 1;
+  const int num_blocks = diy_comm.size();
+  const int magic_k = 8;
+
+  vtkmdiy::Master master(diy_comm, num_threads, -1, 0, [](void* b) {
+    ImageBlock<ImageType>* block = reinterpret_cast<ImageBlock<ImageType>*>(b);
+    delete block;
+  });
+
+  // create an assigner with one block per rank
+  vtkmdiy::ContiguousAssigner assigner(num_blocks, num_blocks);
+  vtkm::rendering::compositing::AddImageBlock<ImageType> create(master, image);
+  const int num_dims = 2;
+  vtkmdiy::RegularDecomposer<vtkmdiy::DiscreteBounds> decomposer(
+    num_dims, global_bounds, num_blocks);
+  decomposer.decompose(diy_comm.rank(), assigner, create);
+  vtkmdiy::RegularSwapPartners partners(decomposer, magic_k,
+                                        false); // false == distance halving
+  vtkmdiy::reduce(master, assigner, partners, reduce_images<ImageType>);
+
+
+  //MPICollect(image, diy_comm);
+  vtkmdiy::all_to_all(
+    master, assigner, vtkm::rendering::compositing::CollectImages<ImageType>(decomposer), magic_k);
+
+  if (diy_comm.rank() == 0)
+  {
+    master.prof.output(m_timing_log);
+  }
+}
+
+void RadixKCompositor::CompositeSurface(vtkmdiy::mpi::communicator& diy_comm,
+                                        vtkm::rendering::compositing::Image& image)
+{
+  CompositeImpl(diy_comm, image);
+}
+
+void RadixKCompositor::CompositeSurface(vtkmdiy::mpi::communicator& diy_comm,
+                                        vtkm::rendering::compositing::PayloadImage& image)
+{
+  CompositeImpl(diy_comm, image);
+}
+
+std::string RadixKCompositor::GetTimingString()
+{
+  std::string res(m_timing_log.str());
+  m_timing_log.str("");
+  return res;
+}
+
+}
+}
+} //namespace vtkm::rendering::compositing

vtkm/rendering/compositing/RadixKCompositor.h

@@ -0,0 +1,54 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_rendering_compositing_RadixKCompositor_h
+#define vtk_m_rendering_compositing_RadixKCompositor_h
+
+#include <vtkm/rendering/vtkm_rendering_export.h>
+
+#include <vtkm/rendering/compositing/Image.h>
+#include <vtkm/rendering/compositing/PayloadImage.h>
+
+#include <vtkm/thirdparty/diy/diy.h>
+#ifdef VTKM_ENABLE_MPI
+//#include <mpi.h>
+//#include <vtkm/thirdparty/diy/mpi.h>
+//#include <vtkm/thirdparty/diy/mpi-cast.h>
+#endif
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+class VTKM_RENDERING_EXPORT RadixKCompositor
+{
+public:
+  RadixKCompositor();
+  ~RadixKCompositor();
+  void CompositeSurface(vtkmdiy::mpi::communicator& diy_comm, Image& image);
+  void CompositeSurface(vtkmdiy::mpi::communicator& diy_comm, PayloadImage& image);
+
+  template <typename ImageType>
+  void CompositeImpl(vtkmdiy::mpi::communicator& diy_comm, ImageType& image);
+
+  std::string GetTimingString();
+
+private:
+  std::stringstream m_timing_log;
+};
+
+}
+}
+} //namespace vtkm::rendering::compositing
+
+#endif //vtk_m_rendering_compositing_RadixKCompositor_h

vtkm/rendering/compositing/VolumePartial.h

@@ -0,0 +1,99 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtkm_rendering_compositing_VolumePartial_h
+#define vtkm_rendering_compositing_VolumePartial_h
+
+#include <limits>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+template <typename FloatType>
+struct VolumePartial
+{
+  typedef FloatType ValueType;
+  int m_pixel_id;
+  float m_depth;
+  float m_pixel[3];
+  float m_alpha;
+
+  VolumePartial()
+    : m_pixel_id(0)
+    , m_depth(0.f)
+    , m_alpha(0.f)
+  {
+    m_pixel[0] = 0;
+    m_pixel[1] = 0;
+    m_pixel[2] = 0;
+  }
+
+  void print() const
+  {
+    std::cout << "[id : " << m_pixel_id << ", red : " << m_pixel[0] << ","
+              << " green : " << m_pixel[1] << ", blue : " << m_pixel[2] << ", alpha " << m_alpha
+              << ", depth : " << m_depth << "]\n";
+  }
+
+  bool operator<(const VolumePartial& other) const
+  {
+    if (m_pixel_id != other.m_pixel_id)
+    {
+      return m_pixel_id < other.m_pixel_id;
+    }
+    else
+    {
+      return m_depth < other.m_depth;
+    }
+  }
+
+  inline void blend(const VolumePartial& other)
+  {
+    if (m_alpha >= 1.f || other.m_alpha == 0.f)
+      return;
+    const float opacity = (1.f - m_alpha);
+    m_pixel[0] += opacity * other.m_pixel[0];
+    m_pixel[1] += opacity * other.m_pixel[1];
+    m_pixel[2] += opacity * other.m_pixel[2];
+    m_alpha += opacity * other.m_alpha;
+    m_alpha = m_alpha > 1.f ? 1.f : m_alpha;
+  }
+
+  static void composite_background(std::vector<VolumePartial>& partials,
+                                   const std::vector<FloatType>& background)
+  {
+    VolumePartial bg_color;
+    bg_color.m_pixel[0] = static_cast<float>(background[0]);
+    bg_color.m_pixel[1] = static_cast<float>(background[1]);
+    bg_color.m_pixel[2] = static_cast<float>(background[2]);
+    bg_color.m_alpha = static_cast<float>(background[3]);
+    //
+    // Gather the unique pixels into the output
+    //
+    const int total_pixels = static_cast<int>(partials.size());
+#ifdef VTKH_OPENMP_ENABLED
+#pragma omp parallel for
+#endif
+    for (int i = 0; i < total_pixels; ++i)
+    {
+      partials[i].blend(bg_color);
+    }
+  }
+};
+
+}
+}
+} //vtkm::rendering::compositing
+
+#endif //vtkm_rendering_compositing_VolumePartial_h

vtkm/rendering/compositing/vtkm_diy_collect.h

@@ -0,0 +1,88 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_rendering_compositing_vtkm_diy_collect_h
+#define vtk_m_rendering_compositing_vtkm_diy_collect_h
+
+#include <vtkm/rendering/vtkm_rendering_export.h>
+
+#include <vtkm/rendering/compositing/Image.h>
+#include <vtkm/rendering/compositing/vtkm_diy_image_block.h>
+#include <vtkmdiy/master.hpp>
+#include <vtkmdiy/partners/swap.hpp>
+#include <vtkmdiy/reduce-operations.hpp>
+#include <vtkmdiy/reduce.hpp>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+template <typename ImageType>
+struct CollectImages
+{
+  const vtkmdiy::RegularDecomposer<vtkmdiy::DiscreteBounds>& m_decomposer;
+
+  CollectImages(const vtkmdiy::RegularDecomposer<vtkmdiy::DiscreteBounds>& decomposer)
+    : m_decomposer(decomposer)
+  {
+  }
+
+  void operator()(void* b, const vtkmdiy::ReduceProxy& proxy) const
+  {
+    ImageBlock<ImageType>* block = reinterpret_cast<ImageBlock<ImageType>*>(b);
+    //
+    // first round we have no incoming. Take the images we have
+    // and sent them to to the right rank
+    //
+    const int collection_rank = 0;
+    if (proxy.in_link().size() == 0)
+    {
+
+      if (proxy.gid() != collection_rank)
+      {
+        int dest_gid = collection_rank;
+        vtkmdiy::BlockID dest = proxy.out_link().target(dest_gid);
+
+        proxy.enqueue(dest, block->m_image);
+        block->m_image.Clear();
+      }
+    } // if
+    else if (proxy.gid() == collection_rank)
+    {
+      ImageType final_image;
+      final_image.InitOriginal(block->m_image);
+      block->m_image.SubsetTo(final_image);
+
+      for (int i = 0; i < proxy.in_link().size(); ++i)
+      {
+        int gid = proxy.in_link().target(i).gid;
+
+        if (gid == collection_rank)
+        {
+          continue;
+        }
+        ImageType incoming;
+        proxy.dequeue(gid, incoming);
+        incoming.SubsetTo(final_image);
+      } // for
+      block->m_image.Swap(final_image);
+    } // else
+
+  } // operator
+};
+
+}
+}
+} //namespace vtkm::rendering::compositing
+
+#endif //vtk_m_rendering_compositing_vtkm_diy_collect_h

vtkm/rendering/compositing/vtkm_diy_image_block.h

@@ -0,0 +1,194 @@
+#ifndef VTKH_DIY_IMAGE_BLOCK_HPP
+#define VTKH_DIY_IMAGE_BLOCK_HPP
+
+#include <vtkm/rendering/compositing/Image.h>
+#include <vtkm/rendering/compositing/PayloadImage.h>
+#include <vtkmdiy/master.hpp>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+template <typename ImageType>
+struct ImageBlock
+{
+  ImageType& m_image;
+  ImageBlock(ImageType& image)
+    : m_image(image)
+  {
+  }
+};
+
+struct MultiImageBlock
+{
+  std::vector<vtkm::rendering::compositing::Image>& m_images;
+  vtkm::rendering::compositing::Image& m_output;
+  MultiImageBlock(std::vector<vtkm::rendering::compositing::Image>& images,
+                  vtkm::rendering::compositing::Image& output)
+    : m_images(images)
+    , m_output(output)
+  {
+  }
+};
+
+template <typename ImageType>
+struct AddImageBlock
+{
+  ImageType& m_image;
+  const vtkmdiy::Master& m_master;
+
+  AddImageBlock(vtkmdiy::Master& master, ImageType& image)
+    : m_image(image)
+    , m_master(master)
+  {
+  }
+  template <typename BoundsType, typename LinkType>
+  void operator()(int gid,
+                  const BoundsType&, // local_bounds
+                  const BoundsType&, // local_with_ghost_bounds
+                  const BoundsType&, // domain_bounds
+                  const LinkType& link) const
+  {
+    ImageBlock<ImageType>* block = new ImageBlock<ImageType>(m_image);
+    LinkType* linked = new LinkType(link);
+    vtkmdiy::Master& master = const_cast<vtkmdiy::Master&>(m_master);
+    master.add(gid, block, linked);
+  }
+};
+
+struct AddMultiImageBlock
+{
+  const vtkmdiy::Master& m_master;
+  std::vector<vtkm::rendering::compositing::Image>& m_images;
+  vtkm::rendering::compositing::Image& m_output;
+
+  AddMultiImageBlock(vtkmdiy::Master& master,
+                     std::vector<vtkm::rendering::compositing::Image>& images,
+                     vtkm::rendering::compositing::Image& output)
+    : m_master(master)
+    , m_images(images)
+    , m_output(output)
+  {
+  }
+  template <typename BoundsType, typename LinkType>
+  void operator()(int gid,
+                  const BoundsType&, // local_bounds
+                  const BoundsType&, // local_with_ghost_bounds
+                  const BoundsType&, // domain_bounds
+                  const LinkType& link) const
+  {
+    MultiImageBlock* block = new MultiImageBlock(m_images, m_output);
+    LinkType* linked = new LinkType(link);
+    vtkmdiy::Master& master = const_cast<vtkmdiy::Master&>(m_master);
+    int lid = master.add(gid, block, linked);
+  }
+};
+
+}
+}
+} //namespace vtkm::rendering::compositing
+
+namespace vtkmdiy
+{
+
+template <>
+struct Serialization<vtkm::rendering::compositing::PayloadImage>
+{
+  static void save(BinaryBuffer& bb, const vtkm::rendering::compositing::PayloadImage& image)
+  {
+    vtkmdiy::save(bb, image.OrigBounds.X.Min);
+    vtkmdiy::save(bb, image.OrigBounds.Y.Min);
+    vtkmdiy::save(bb, image.OrigBounds.Z.Min);
+    vtkmdiy::save(bb, image.OrigBounds.X.Max);
+    vtkmdiy::save(bb, image.OrigBounds.Y.Max);
+    vtkmdiy::save(bb, image.OrigBounds.Z.Max);
+
+    vtkmdiy::save(bb, image.Bounds.X.Min);
+    vtkmdiy::save(bb, image.Bounds.Y.Min);
+    vtkmdiy::save(bb, image.Bounds.Z.Min);
+    vtkmdiy::save(bb, image.Bounds.X.Max);
+    vtkmdiy::save(bb, image.Bounds.Y.Max);
+    vtkmdiy::save(bb, image.Bounds.Z.Max);
+
+    vtkmdiy::save(bb, image.Payloads);
+    vtkmdiy::save(bb, image.PayloadBytes);
+    vtkmdiy::save(bb, image.Depths);
+    vtkmdiy::save(bb, image.OrigRank);
+  }
+
+  static void load(BinaryBuffer& bb, vtkm::rendering::compositing::PayloadImage& image)
+  {
+    vtkmdiy::load(bb, image.OrigBounds.X.Min);
+    vtkmdiy::load(bb, image.OrigBounds.Y.Min);
+    vtkmdiy::load(bb, image.OrigBounds.Z.Min);
+    vtkmdiy::load(bb, image.OrigBounds.X.Max);
+    vtkmdiy::load(bb, image.OrigBounds.Y.Max);
+    vtkmdiy::load(bb, image.OrigBounds.Z.Max);
+
+    vtkmdiy::load(bb, image.Bounds.X.Min);
+    vtkmdiy::load(bb, image.Bounds.Y.Min);
+    vtkmdiy::load(bb, image.Bounds.Z.Min);
+    vtkmdiy::load(bb, image.Bounds.X.Max);
+    vtkmdiy::load(bb, image.Bounds.Y.Max);
+    vtkmdiy::load(bb, image.Bounds.Z.Max);
+
+    vtkmdiy::load(bb, image.Payloads);
+    vtkmdiy::load(bb, image.PayloadBytes);
+    vtkmdiy::load(bb, image.Depths);
+    vtkmdiy::load(bb, image.OrigRank);
+  }
+};
+
+template <>
+struct Serialization<vtkm::rendering::compositing::Image>
+{
+  static void save(BinaryBuffer& bb, const vtkm::rendering::compositing::Image& image)
+  {
+    vtkmdiy::save(bb, image.OrigBounds.X.Min);
+    vtkmdiy::save(bb, image.OrigBounds.Y.Min);
+    vtkmdiy::save(bb, image.OrigBounds.Z.Min);
+    vtkmdiy::save(bb, image.OrigBounds.X.Max);
+    vtkmdiy::save(bb, image.OrigBounds.Y.Max);
+    vtkmdiy::save(bb, image.OrigBounds.Z.Max);
+
+    vtkmdiy::save(bb, image.Bounds.X.Min);
+    vtkmdiy::save(bb, image.Bounds.Y.Min);
+    vtkmdiy::save(bb, image.Bounds.Z.Min);
+    vtkmdiy::save(bb, image.Bounds.X.Max);
+    vtkmdiy::save(bb, image.Bounds.Y.Max);
+    vtkmdiy::save(bb, image.Bounds.Z.Max);
+
+    vtkmdiy::save(bb, image.Pixels);
+    vtkmdiy::save(bb, image.Depths);
+    vtkmdiy::save(bb, image.OrigRank);
+    vtkmdiy::save(bb, image.CompositeOrder);
+  }
+
+  static void load(BinaryBuffer& bb, vtkm::rendering::compositing::Image& image)
+  {
+    vtkmdiy::load(bb, image.OrigBounds.X.Min);
+    vtkmdiy::load(bb, image.OrigBounds.Y.Min);
+    vtkmdiy::load(bb, image.OrigBounds.Z.Min);
+    vtkmdiy::load(bb, image.OrigBounds.X.Max);
+    vtkmdiy::load(bb, image.OrigBounds.Y.Max);
+    vtkmdiy::load(bb, image.OrigBounds.Z.Max);
+
+    vtkmdiy::load(bb, image.Bounds.X.Min);
+    vtkmdiy::load(bb, image.Bounds.Y.Min);
+    vtkmdiy::load(bb, image.Bounds.Z.Min);
+    vtkmdiy::load(bb, image.Bounds.X.Max);
+    vtkmdiy::load(bb, image.Bounds.Y.Max);
+    vtkmdiy::load(bb, image.Bounds.Z.Max);
+
+    vtkmdiy::load(bb, image.Pixels);
+    vtkmdiy::load(bb, image.Depths);
+    vtkmdiy::load(bb, image.OrigRank);
+    vtkmdiy::load(bb, image.CompositeOrder);
+  }
+};
+} //namespace vtkmdiy
+
+#endif

vtkm/rendering/compositing/vtkm_diy_partial_blocks.h

@@ -0,0 +1,198 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtkm_rendering_compositing_vtkm_diy_partial_blocks_h
+#define vtkm_rendering_compositing_vtkm_diy_partial_blocks_h
+
+#include <vtkm/thirdparty/diy/master.h>
+
+#include <vtkm/rendering/compositing/AbsorptionPartial.h>
+#include <vtkm/rendering/compositing/EmissionPartial.h>
+#include <vtkm/rendering/compositing/VolumePartial.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+//--------------------------------------Volume Block Structure-----------------------------------
+template <typename FloatType>
+struct VolumeBlock
+{
+  typedef vtkmdiy::DiscreteBounds Bounds;
+  typedef VolumePartial<FloatType> PartialType;
+  std::vector<VolumePartial<FloatType>>& m_partials;
+  VolumeBlock(std::vector<VolumePartial<FloatType>>& partials)
+    : m_partials(partials)
+  {
+  }
+};
+
+
+//--------------------------------------Absorption Block Structure------------------------------
+template <typename FloatType>
+struct AbsorptionBlock
+{
+  typedef vtkmdiy::DiscreteBounds Bounds;
+  typedef AbsorptionPartial<FloatType> PartialType;
+  std::vector<AbsorptionPartial<FloatType>>& m_partials;
+
+  AbsorptionBlock(std::vector<AbsorptionPartial<FloatType>>& partials)
+    : m_partials(partials)
+  {
+  }
+};
+
+//--------------------------------------Emission Block Structure------------------------------
+template <typename FloatType>
+struct EmissionBlock
+{
+  typedef vtkmdiy::DiscreteBounds Bounds;
+  typedef EmissionPartial<FloatType> PartialType;
+  std::vector<EmissionPartial<FloatType>>& m_partials;
+
+  EmissionBlock(std::vector<EmissionPartial<FloatType>>& partials)
+    : m_partials(partials)
+  {
+  }
+};
+
+//--------------------------------------Add Block Template-----------------------------------
+template <typename BlockType>
+struct AddBlock
+{
+  typedef typename BlockType::PartialType PartialType;
+  typedef BlockType Block;
+  std::vector<PartialType>& m_partials;
+  const vtkmdiy::Master& m_master;
+
+  AddBlock(vtkmdiy::Master& master, std::vector<PartialType>& partials)
+    : m_master(master)
+    , m_partials(partials)
+  {
+  }
+  template <typename BoundsType, typename LinkType>
+  void operator()(int gid,
+                  const BoundsType& local_bounds,
+                  const BoundsType& local_with_ghost_bounds,
+                  const BoundsType& domain_bounds,
+                  const LinkType& link) const
+  {
+    (void)local_bounds;
+    (void)domain_bounds;
+    (void)local_with_ghost_bounds;
+    Block* block = new Block(m_partials);
+    LinkType* rg_link = new LinkType(link);
+    vtkmdiy::Master& master = const_cast<vtkmdiy::Master&>(m_master);
+    int lid = master.add(gid, block, rg_link);
+    (void)lid;
+  }
+};
+
+}
+}
+} //vtkm::rendering::compositing
+
+//-------------------------------Serialization Specializations--------------------------------
+namespace vtkmdiy
+{
+
+template <>
+struct Serialization<vtkm::rendering::compositing::AbsorptionPartial<double>>
+{
+
+  static void save(BinaryBuffer& bb,
+                   const vtkm::rendering::compositing::AbsorptionPartial<double>& partial)
+  {
+    vtkmdiy::save(bb, partial.m_bins);
+    vtkmdiy::save(bb, partial.m_pixel_id);
+    vtkmdiy::save(bb, partial.m_depth);
+  }
+
+  static void load(vtkmdiy::BinaryBuffer& bb,
+                   vtkm::rendering::compositing::AbsorptionPartial<double>& partial)
+  {
+    vtkmdiy::load(bb, partial.m_bins);
+    vtkmdiy::load(bb, partial.m_pixel_id);
+    vtkmdiy::load(bb, partial.m_depth);
+  }
+};
+
+template <>
+struct Serialization<vtkm::rendering::compositing::AbsorptionPartial<float>>
+{
+
+  static void save(BinaryBuffer& bb,
+                   const vtkm::rendering::compositing::AbsorptionPartial<float>& partial)
+  {
+    vtkmdiy::save(bb, partial.m_bins);
+    vtkmdiy::save(bb, partial.m_pixel_id);
+    vtkmdiy::save(bb, partial.m_depth);
+  }
+
+  static void load(BinaryBuffer& bb,
+                   vtkm::rendering::compositing::AbsorptionPartial<float>& partial)
+  {
+    vtkmdiy::load(bb, partial.m_bins);
+    vtkmdiy::load(bb, partial.m_pixel_id);
+    vtkmdiy::load(bb, partial.m_depth);
+  }
+};
+
+template <>
+struct Serialization<vtkm::rendering::compositing::EmissionPartial<double>>
+{
+
+  static void save(BinaryBuffer& bb,
+                   const vtkm::rendering::compositing::EmissionPartial<double>& partial)
+  {
+    vtkmdiy::save(bb, partial.m_bins);
+    vtkmdiy::save(bb, partial.m_emission_bins);
+    vtkmdiy::save(bb, partial.m_pixel_id);
+    vtkmdiy::save(bb, partial.m_depth);
+  }
+
+  static void load(BinaryBuffer& bb, vtkm::rendering::compositing::EmissionPartial<double>& partial)
+  {
+    vtkmdiy::load(bb, partial.m_bins);
+    vtkmdiy::load(bb, partial.m_emission_bins);
+    vtkmdiy::load(bb, partial.m_pixel_id);
+    vtkmdiy::load(bb, partial.m_depth);
+  }
+};
+
+template <>
+struct Serialization<vtkm::rendering::compositing::EmissionPartial<float>>
+{
+
+  static void save(BinaryBuffer& bb,
+                   const vtkm::rendering::compositing::EmissionPartial<float>& partial)
+  {
+    vtkmdiy::save(bb, partial.m_bins);
+    vtkmdiy::save(bb, partial.m_emission_bins);
+    vtkmdiy::save(bb, partial.m_pixel_id);
+    vtkmdiy::save(bb, partial.m_depth);
+  }
+
+  static void load(BinaryBuffer& bb, vtkm::rendering::compositing::EmissionPartial<float>& partial)
+  {
+    vtkmdiy::load(bb, partial.m_bins);
+    vtkmdiy::load(bb, partial.m_emission_bins);
+    vtkmdiy::load(bb, partial.m_pixel_id);
+    vtkmdiy::load(bb, partial.m_depth);
+  }
+};
+
+} //vtkmdiy
+
+
+#endif //vtkm_rendering_compositing_vtkm_diy_partial_blocks_h

vtkm/rendering/compositing/vtkm_diy_partial_collect.h

@@ -0,0 +1,215 @@
+//============================================================================
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtkm_rendering_compositing_vtkm_diy_partial_collect_h
+#define vtkm_rendering_compositing_vtkm_diy_partial_collect_h
+
+#include <vtkm/rendering/compositing/AbsorptionPartial.h>
+#include <vtkm/rendering/compositing/EmissionPartial.h>
+#include <vtkm/rendering/compositing/VolumePartial.h>
+#include <vtkm/rendering/compositing/vtkm_diy_partial_blocks.h>
+
+#include <vtkm/thirdparty/diy/assigner.h>
+#include <vtkm/thirdparty/diy/decomposition.h>
+#include <vtkm/thirdparty/diy/diy.h>
+#include <vtkm/thirdparty/diy/master.h>
+#include <vtkm/thirdparty/diy/mpi-cast.h>
+#include <vtkm/thirdparty/diy/reduce-operations.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+//
+// Collect struct sends all data to a single node.
+//
+template <typename BlockType>
+struct Collect
+{
+  const vtkmdiy::RegularDecomposer<vtkmdiy::ContinuousBounds>& m_decomposer;
+
+  Collect(const vtkmdiy::RegularDecomposer<vtkmdiy::ContinuousBounds>& decomposer)
+    : m_decomposer(decomposer)
+  {
+  }
+
+  void operator()(void* v_block, const vtkmdiy::ReduceProxy& proxy) const
+  {
+    BlockType* block = static_cast<BlockType*>(v_block);
+    //
+    // first round we have no incoming. Take the partials we have
+    // and sent them to to the right rank
+    //
+    const int collection_rank = 0;
+    if (proxy.in_link().size() == 0 && proxy.gid() != collection_rank)
+    {
+      int dest_gid = collection_rank;
+      vtkmdiy::BlockID dest = proxy.out_link().target(dest_gid);
+      proxy.enqueue(dest, block->m_partials);
+
+      block->m_partials.clear();
+
+    } // if
+    else if (proxy.gid() == collection_rank)
+    {
+
+      for (int i = 0; i < proxy.in_link().size(); ++i)
+      {
+        int gid = proxy.in_link().target(i).gid;
+        if (gid == collection_rank)
+        {
+          continue;
+        }
+        //TODO: leave the paritals that start here, here
+        std::vector<typename BlockType::PartialType> incoming_partials;
+        proxy.dequeue(gid, incoming_partials);
+        const int incoming_size = incoming_partials.size();
+        // TODO: make this a std::copy
+        for (int j = 0; j < incoming_size; ++j)
+        {
+          block->m_partials.push_back(incoming_partials[j]);
+        }
+      } // for
+    }   // else
+
+  } // operator
+};
+
+//
+// collect uses the all-to-all construct to perform a gather to
+// the root rank. All other ranks will have no data
+//
+template <typename AddBlockType>
+void collect_detail(std::vector<typename AddBlockType::PartialType>& partials, MPI_Comm comm)
+{
+  typedef typename AddBlockType::Block Block;
+
+  vtkmdiy::mpi::communicator world(vtkmdiy::mpi::make_DIY_MPI_Comm(comm));
+  std::cout << __FILE__ << " " << __LINE__ << std::endl;
+  std::cout << "             DRP: Is this the right dimension???" << std::endl;
+  vtkmdiy::ContinuousBounds global_bounds(1); //DRP???
+  global_bounds.min[0] = 0;
+  global_bounds.max[0] = 1;
+
+  // tells diy to use all availible threads
+  const int num_threads = -1;
+  const int num_blocks = world.size();
+  const int magic_k = 2;
+
+  vtkmdiy::Master master(world, num_threads);
+
+  // create an assigner with one block per rank
+  vtkmdiy::ContiguousAssigner assigner(num_blocks, num_blocks);
+  AddBlockType create(master, partials);
+
+  const int dims = 1;
+  vtkmdiy::RegularDecomposer<vtkmdiy::ContinuousBounds> decomposer(dims, global_bounds, num_blocks);
+  decomposer.decompose(world.rank(), assigner, create);
+
+  vtkmdiy::all_to_all(master, assigner, Collect<Block>(decomposer), magic_k);
+}
+
+template <typename T>
+void collect(std::vector<T>& partials, MPI_Comm comm);
+
+template <>
+void collect<VolumePartial<float>>(std::vector<VolumePartial<float>>& partials, MPI_Comm comm)
+{
+  collect_detail<vtkm::rendering::compositing::AddBlock<VolumeBlock<float>>>(partials, comm);
+}
+
+
+template <>
+void collect<vtkm::rendering::compositing::VolumePartial<double>>(
+  std::vector<vtkm::rendering::compositing::VolumePartial<double>>& partials,
+  MPI_Comm comm)
+{
+  collect_detail<
+    vtkm::rendering::compositing::AddBlock<vtkm::rendering::compositing::VolumeBlock<double>>>(
+    partials, comm);
+}
+
+template <>
+void collect<vtkm::rendering::compositing::AbsorptionPartial<double>>(
+  std::vector<vtkm::rendering::compositing::AbsorptionPartial<double>>& partials,
+  MPI_Comm comm)
+{
+  collect_detail<
+    vtkm::rendering::compositing::AddBlock<vtkm::rendering::compositing::AbsorptionBlock<double>>>(
+    partials, comm);
+}
+
+template <>
+void collect<vtkm::rendering::compositing::AbsorptionPartial<float>>(
+  std::vector<vtkm::rendering::compositing::AbsorptionPartial<float>>& partials,
+  MPI_Comm comm)
+{
+  collect_detail<
+    vtkm::rendering::compositing::AddBlock<vtkm::rendering::compositing::AbsorptionBlock<float>>>(
+    partials, comm);
+}
+
+template <>
+void collect<vtkm::rendering::compositing::EmissionPartial<double>>(
+  std::vector<vtkm::rendering::compositing::EmissionPartial<double>>& partials,
+  MPI_Comm comm)
+{
+  collect_detail<vtkm::rendering::compositing::AddBlock<EmissionBlock<double>>>(partials, comm);
+}
+
+template <>
+void collect<vtkm::rendering::compositing::EmissionPartial<float>>(
+  std::vector<vtkm::rendering::compositing::EmissionPartial<float>>& partials,
+  MPI_Comm comm)
+{
+  collect_detail<
+    vtkm::rendering::compositing::AddBlock<vtkm::rendering::compositing::EmissionBlock<float>>>(
+    partials, comm);
+}
+
+}
+}
+} //vtkm::rendering::compositing
+
+#endif //vtkm_rendering_compositing_vtkm_diy_partial_collect_h
+
+
+#if 0
+  /*
+template <>
+void collect<vtkm::rendering::compositing::EmissionPartial<float>>(std::vector<vtkm::rendering::compositing::EmissionPartial<float>>& partials, MPI_Comm comm)
+{
+  collect_detail<vtkm::rendering::compositing::AddBlock<vtkm::rendering::compositing::EmissionBlock<float>>>(partials, comm);
+}
+  */
+  }
+  * /
+
+  }
+
+
+}
+
+/*
+
+#endif //vtkm_rendering_compositing_vtkm_diy_partial_collect_h
+template <>
+template <>
+void colloct<EmissioiPartial<float>>(std::vector<Em ssionPartial<cloat>>& partials, MPI_Comm comm)
+{
+  ocollect_detail < ect<::missionPa::tial<float>::AddBloct<EdissvonBlock<float>>>(tor < Emis, simm);
+}al<float>>& partials, MPI_Comm comm)
+{
+  collect_detail<vtkm::rendering::compositing::AddBlock<EmissionBlock<float>>>(partials, comm);
+}
+
+*/
+#endif

vtkm/rendering/compositing/vtkm_diy_partial_redistribute.h

@@ -0,0 +1,205 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtkm_rendering_compositing_vtkm_diy_partial_redistribute_h
+#define vtkm_rendering_compositing_vtkm_diy_partial_redistribute_h
+
+#include <map>
+#include <vtkm/rendering/compositing/vtkm_diy_partial_blocks.h>
+#include <vtkm/thirdparty/diy/assigner.h>
+#include <vtkm/thirdparty/diy/decomposition.h>
+#include <vtkm/thirdparty/diy/master.h>
+#include <vtkm/thirdparty/diy/point.h>
+#include <vtkm/thirdparty/diy/reduce-operations.h>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+//
+// Redistributes partial composites to the ranks that owns
+// that sectoon of the image. Currently, the domain is decomposed
+// in 1-D from min_pixel to max_pixel.
+//
+template <typename BlockType>
+struct Redistribute
+{
+  const vtkmdiy::RegularDecomposer<vtkmdiy::DiscreteBounds>& m_decomposer;
+
+  Redistribute(const vtkmdiy::RegularDecomposer<vtkmdiy::DiscreteBounds>& decomposer)
+    : m_decomposer(decomposer)
+  {
+  }
+
+  void operator()(void* v_block, const vtkmdiy::ReduceProxy& proxy) const
+  {
+    BlockType* block = static_cast<BlockType*>(v_block);
+    //
+    // first round we have no incoming. Take the partials we have
+    // and sent them to to the right rank
+    //
+    if (proxy.in_link().size() == 0)
+    {
+      const int size = block->m_partials.size();
+      std::map<vtkmdiy::BlockID, std::vector<typename BlockType::PartialType>> outgoing;
+
+      for (int i = 0; i < size; ++i)
+      {
+        vtkmdiy::Point<int, VTKMDIY_MAX_DIM> point;
+        point[0] = block->m_partials[i].m_pixel_id;
+        int dest_gid = m_decomposer.point_to_gid(point);
+        vtkmdiy::BlockID dest = proxy.out_link().target(dest_gid);
+        outgoing[dest].push_back(block->m_partials[i]);
+      } //for
+
+      block->m_partials.clear();
+
+
+      for (int i = 0; i < proxy.out_link().size(); ++i)
+      {
+        int dest_gid = proxy.out_link().target(i).gid;
+        vtkmdiy::BlockID dest = proxy.out_link().target(dest_gid);
+        proxy.enqueue(dest, outgoing[dest]);
+        //outgoing[dest].clear();
+      }
+
+    } // if
+    else
+    {
+      for (int i = 0; i < proxy.in_link().size(); ++i)
+      {
+        int gid = proxy.in_link().target(i).gid;
+        std::vector<typename BlockType::PartialType> incoming_partials;
+        proxy.dequeue(gid, incoming_partials);
+        const int incoming_size = incoming_partials.size();
+        // TODO: make this a std::copy
+        for (int j = 0; j < incoming_size; ++j)
+        {
+          block->m_partials.push_back(incoming_partials[j]);
+        }
+      } // for
+
+    }                            // else
+    MPI_Barrier(MPI_COMM_WORLD); //HACK
+  }                              // operator
+};
+
+
+template <typename AddBlockType>
+void redistribute_detail(std::vector<typename AddBlockType::PartialType>& partials,
+                         MPI_Comm comm,
+                         const int& domain_min_pixel,
+                         const int& domain_max_pixel)
+{
+  typedef typename AddBlockType::Block Block;
+
+  vtkmdiy::mpi::communicator world(vtkmdiy::mpi::make_DIY_MPI_Comm(comm));
+  std::cout << __FILE__ << " " << __LINE__ << std::endl;
+  std::cout << "             DRP: Is this the right dimension???" << std::endl;
+  vtkmdiy::DiscreteBounds global_bounds(1); //DRP???
+  global_bounds.min[0] = domain_min_pixel;
+  global_bounds.max[0] = domain_max_pixel;
+
+  // tells diy to use all availible threads
+  const int num_threads = 1;
+  const int num_blocks = world.size();
+  const int magic_k = 2;
+
+  vtkmdiy::Master master(world, num_threads);
+
+  // create an assigner with one block per rank
+  vtkmdiy::ContiguousAssigner assigner(num_blocks, num_blocks);
+  AddBlockType create(master, partials);
+
+  const int dims = 1;
+  vtkmdiy::RegularDecomposer<vtkmdiy::DiscreteBounds> decomposer(dims, global_bounds, num_blocks);
+  decomposer.decompose(world.rank(), assigner, create);
+  vtkmdiy::all_to_all(master, assigner, Redistribute<Block>(decomposer), magic_k);
+}
+
+//
+// Define a default template that cannot be instantiated
+//
+template <typename T>
+void redistribute(std::vector<T>& partials,
+                  MPI_Comm comm,
+                  const int& domain_min_pixel,
+                  const int& domain_max_pixel);
+// ----------------------------- VolumePartial Specialization------------------------------------------
+template <>
+void redistribute<VolumePartial<float>>(std::vector<VolumePartial<float>>& partials,
+                                        MPI_Comm comm,
+                                        const int& domain_min_pixel,
+                                        const int& domain_max_pixel)
+{
+  redistribute_detail<AddBlock<VolumeBlock<float>>>(
+    partials, comm, domain_min_pixel, domain_max_pixel);
+}
+
+template <>
+void redistribute<VolumePartial<double>>(std::vector<VolumePartial<double>>& partials,
+                                         MPI_Comm comm,
+                                         const int& domain_min_pixel,
+                                         const int& domain_max_pixel)
+{
+  redistribute_detail<AddBlock<VolumeBlock<double>>>(
+    partials, comm, domain_min_pixel, domain_max_pixel);
+}
+
+// ----------------------------- AbsorpPartial Specialization------------------------------------------
+template <>
+void redistribute<AbsorptionPartial<double>>(std::vector<AbsorptionPartial<double>>& partials,
+                                             MPI_Comm comm,
+                                             const int& domain_min_pixel,
+                                             const int& domain_max_pixel)
+{
+  redistribute_detail<AddBlock<AbsorptionBlock<double>>>(
+    partials, comm, domain_min_pixel, domain_max_pixel);
+}
+
+template <>
+void redistribute<AbsorptionPartial<float>>(std::vector<AbsorptionPartial<float>>& partials,
+                                            MPI_Comm comm,
+                                            const int& domain_min_pixel,
+                                            const int& domain_max_pixel)
+{
+  redistribute_detail<AddBlock<AbsorptionBlock<float>>>(
+    partials, comm, domain_min_pixel, domain_max_pixel);
+}
+
+// ----------------------------- EmissPartial Specialization------------------------------------------
+template <>
+void redistribute<EmissionPartial<double>>(std::vector<EmissionPartial<double>>& partials,
+                                           MPI_Comm comm,
+                                           const int& domain_min_pixel,
+                                           const int& domain_max_pixel)
+{
+  redistribute_detail<AddBlock<EmissionBlock<double>>>(
+    partials, comm, domain_min_pixel, domain_max_pixel);
+}
+
+template <>
+void redistribute<EmissionPartial<float>>(std::vector<EmissionPartial<float>>& partials,
+                                          MPI_Comm comm,
+                                          const int& domain_min_pixel,
+                                          const int& domain_max_pixel)
+{
+  redistribute_detail<AddBlock<EmissionBlock<float>>>(
+    partials, comm, domain_min_pixel, domain_max_pixel);
+}
+
+}
+}
+} //vtkm::rendering::compositing
+
+#endif //vtkm_rendering_compositing_vtkm_diy_partial_redistribute_h

vtkm/rendering/compositing/vtkm_diy_utils.h

@@ -0,0 +1,65 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtkm_rendering_compositing_vtkm_diy_utils_h
+#define vtkm_rendering_compositing_vtkm_diy_utils_h
+
+#include <vtkm/Bounds.h>
+#include <vtkmdiy/decomposition.hpp>
+
+namespace vtkm
+{
+namespace rendering
+{
+namespace compositing
+{
+
+static vtkm::Bounds DIYBoundsToVTKM(const vtkmdiy::DiscreteBounds& bounds)
+{
+  vtkm::Bounds vtkm_bounds;
+
+  vtkm_bounds.X.Min = bounds.min[0];
+  vtkm_bounds.Y.Min = bounds.min[1];
+  vtkm_bounds.Z.Min = bounds.min[2];
+
+  vtkm_bounds.X.Max = bounds.max[0];
+  vtkm_bounds.Y.Max = bounds.max[1];
+  vtkm_bounds.Z.Max = bounds.max[2];
+  return vtkm_bounds;
+}
+
+static vtkmdiy::DiscreteBounds VTKMBoundsToDIY(const vtkm::Bounds& bounds)
+{
+  vtkmdiy::DiscreteBounds diy_bounds(3);
+
+  diy_bounds.min[0] = bounds.X.Min;
+  diy_bounds.min[1] = bounds.Y.Min;
+
+  diy_bounds.max[0] = bounds.X.Max;
+  diy_bounds.max[1] = bounds.Y.Max;
+
+  if (bounds.Z.IsNonEmpty())
+  {
+    diy_bounds.min[2] = bounds.Z.Min;
+    diy_bounds.max[2] = bounds.Z.Max;
+  }
+  else
+  {
+    diy_bounds.min[2] = 0;
+    diy_bounds.max[2] = 0;
+  }
+  return diy_bounds;
+}
+
+}
+}
+} //vtkm::rendering::compositing
+
+#endif //vtkm_rendering_compositing_vtkm_diy_utils_h

vtkm/rendering/raytracing/RayTracer.cxx

@@ -56,16 +56,18 @@ public:
       , LookAt(lookAt)
     {
       //Set up some default lighting parameters for now
-      LightAbmient[0] = .5f;
-      LightAbmient[1] = .5f;
-      LightAbmient[2] = .5f;
-      LightDiffuse[0] = .7f;
-      LightDiffuse[1] = .7f;
-      LightDiffuse[2] = .7f;
-      LightSpecular[0] = .7f;
-      LightSpecular[1] = .7f;
-      LightSpecular[2] = .7f;
-      SpecularExponent = 20.f;
+      LightAbmient[0] = 0.1f;
+      LightAbmient[1] = 0.1f;
+      LightAbmient[2] = 0.1f;
+
+      LightDiffuse[0] = 1.0f;
+      LightDiffuse[1] = 1.0f;
+      LightDiffuse[2] = 1.0f;
+
+      LightSpecular[0] = 1.0f;
+      LightSpecular[1] = 1.0f;
+      LightSpecular[2] = 1.0f;
+      SpecularExponent = 50.0f;
     }
 
     using ControlSignature =
@@ -116,14 +118,30 @@ public:
       // clamp color index
       colorIdx = vtkm::Max(0, colorIdx);
       colorIdx = vtkm::Min(colorMapSize - 1, colorIdx);
-      color = colorMap.Get(colorIdx);
+      vtkm::Vec<Precision, 4> diffuseColor = colorMap.Get(colorIdx);
 
-      color[0] *= vtkm::Min(
-        LightAbmient[0] + LightDiffuse[0] * cosTheta + LightSpecular[0] * specularConstant, one);
-      color[1] *= vtkm::Min(
-        LightAbmient[1] + LightDiffuse[1] * cosTheta + LightSpecular[1] * specularConstant, one);
-      color[2] *= vtkm::Min(
-        LightAbmient[2] + LightDiffuse[2] * cosTheta + LightSpecular[2] * specularConstant, one);
+      // Add ambient lighting
+      color[0] = LightAbmient[0];
+      color[1] = LightAbmient[1];
+      color[2] = LightAbmient[2];
+
+      // Add diffuse lighting
+      color[0] += LightDiffuse[0] * diffuseColor[0] * cosTheta;
+      color[1] += LightDiffuse[1] * diffuseColor[1] * cosTheta;
+      color[2] += LightDiffuse[2] * diffuseColor[2] * cosTheta;
+
+      // Add specular lighting
+      color[0] += LightSpecular[0] * specularConstant;
+      color[1] += LightSpecular[1] * specularConstant;
+      color[2] += LightSpecular[2] * specularConstant;
+
+      // Set alpha to 1.0 (opaque)
+      color[3] = 1.0f;
+
+      // Clamp color values to [0,1]
+      color[0] = vtkm::Clamp(color[0], zero, one);
+      color[1] = vtkm::Clamp(color[1], zero, one);
+      color[2] = vtkm::Clamp(color[2], zero, one);
 
       colors.Set(offset + 0, color[0]);
       colors.Set(offset + 1, color[1]);

vtkm/rendering/raytracing/RayTracer.h

@@ -34,6 +34,7 @@ protected:
   vtkm::cont::ArrayHandle<vtkm::Vec4f_32> ColorMap;
   vtkm::Range ScalarRange;
   bool Shade;
+  vtkm::Vec3f_32 LightPosition;
 
   template <typename Precision>
   void RenderOnDevice(Ray<Precision>& rays);

vtkm/rendering/testing/CMakeLists.txt

@@ -25,4 +25,18 @@ set(unit_tests
   UnitTestMapperGlyphVector.cxx
 )
 
-vtkm_unit_tests(SOURCES ${unit_tests})
+vtkm_unit_tests(SOURCES ${unit_tests} LIBRARIES vtkm_source vtkm_filter_field_transform)
+
+#add distributed tests i.e.test to run with MPI
+#if MPI is enabled.
+if (VTKm_ENABLE_MPI)
+  set(mpi_unit_tests
+    UnitTestImageCompositing.cxx
+  )
+  vtkm_unit_tests(
+    MPI
+    DEVICE_SOURCES ${mpi_unit_tests}
+    LIBRARIES vtkm_source
+    USE_VTKM_JOB_POOL
+  )
+endif()

vtkm/rendering/testing/UnitTestImageCompositing.cxx

@@ -0,0 +1,264 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#include <vtkm/cont/testing/Testing.h>
+#include <vtkm/rendering/compositing/Compositor.h>
+#include <vtkm/rendering/compositing/Image.h>
+#include <vtkm/rendering/compositing/ImageCompositor.h>
+#include <vtkm/rendering/testing/t_vtkm_test_utils.h>
+#include <vtkm/source/Tangle.h>
+
+#include <vtkm/rendering/Actor.h>
+#include <vtkm/rendering/CanvasRayTracer.h>
+#include <vtkm/rendering/MapperConnectivity.h>
+#include <vtkm/rendering/MapperRayTracer.h>
+#include <vtkm/rendering/MapperVolume.h>
+#include <vtkm/rendering/MapperWireframer.h>
+#include <vtkm/rendering/Scene.h>
+#include <vtkm/rendering/View3D.h>
+#include <vtkm/rendering/testing/RenderTest.h>
+
+#include <vtkm/io/VTKDataSetReader.h>
+
+#include <vtkm/filter/clean_grid/CleanGrid.h>
+
+namespace
+{
+
+template <typename T>
+T* GetVTKMPointer(vtkm::cont::ArrayHandle<T>& handle)
+{
+  return handle.WritePortal().GetArray();
+}
+
+
+vtkm::cont::DataSet ReadDS(int rank)
+{
+
+  std::string vtkFile;
+
+  vtkm::io::VTKDataSetReader reader(vtkFile);
+}
+
+#if 0
+vtkm::rendering::compositing::Image ConstImage(const std::size_t& width,
+                                               const std::size_t& height,
+                                               const vtkm::Vec4f& rgba,
+                                               const vtkm::FloatDefault& depth)
+{
+  auto numPix = width * height;
+  std::vector<vtkm::FloatDefault> rgbaVals(numPix * 4);
+  std::vector<vtkm::FloatDefault> depthVals(numPix, depth);
+
+  for (std::size_t i = 0; i < numPix; i++)
+  {
+    rgbaVals[i * 4 + 0] = rgba[0];
+    rgbaVals[i * 4 + 1] = rgba[1];
+    rgbaVals[i * 4 + 2] = rgba[2];
+    rgbaVals[i * 4 + 3] = rgba[3];
+  }
+
+  vtkm::rendering::compositing::Image img(vtkm::Bounds(0, width, 0, height, 0, 1));
+  img.Init(rgbaVals.data(), depthVals.data(), width, height);
+
+  return img;
+}
+#endif
+
+void TestImageComposite()
+{
+#if 0
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+
+  std::size_t width = 4, height = 4;
+
+  //res is the background, initially black.
+  auto img0 = ConstImage(width, height, { 1, 0, 0, 1 }, 1.0);
+
+  auto img1 = ConstImage(width, height, { 0, 1, 1, .5 }, 0.5);
+
+  vtkm::rendering::compositing::Compositor compositor;
+
+  compositor.SetCompositeMode(vtkm::rendering::compositing::Compositor::Z_BUFFER_SURFACE);
+  vtkm::rendering::compositing::Image img;
+
+  if (comm.rank() == 0)
+    img = ConstImage(width, height, { 1, 0, 0, 1 }, 1.0);
+  else
+    img = ConstImage(width, height, { 0, 1, 1, .5 }, 0.5);
+
+  compositor.AddImage(img.m_pixels.data(), img.m_depths.data(), width, height);
+
+  auto res = compositor.Composite();
+
+  //vtkm::rendering::compositing::ImageCompositor imgCompositor;
+  //compositor.ZBufferComposite(res, img);
+  //compositor.Blend(res, img);
+
+  if (comm.rank() == 0)
+  {
+    for (int i = 0; i < width * height; i++)
+    {
+      std::cout << i << ": ";
+      std::cout << (int)res.m_pixels[i * 4 + 0] << " ";
+      std::cout << (int)res.m_pixels[i * 4 + 1] << " ";
+      std::cout << (int)res.m_pixels[i * 4 + 2] << " ";
+      std::cout << (int)res.m_pixels[i * 4 + 3] << " ";
+      std::cout << res.m_depths[i] << std::endl;
+    }
+  }
+#endif
+}
+
+void TestRenderComposite()
+{
+  using vtkm::rendering::CanvasRayTracer;
+  using vtkm::rendering::MapperRayTracer;
+  using vtkm::rendering::MapperVolume;
+  using vtkm::rendering::MapperWireframer;
+
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+
+  int numBlocks = comm.size() * 1;
+  int rank = comm.rank();
+  int dsPerRank = 2;
+
+  vtkm::rendering::Camera camera;
+  camera.SetLookAt(vtkm::Vec3f_32(1.0, 0.5, 0.5));
+  camera.SetViewUp(vtkm::make_Vec(0.f, 1.f, 0.f));
+  camera.SetClippingRange(1.f, 10.f);
+  camera.SetFieldOfView(60.f);
+  camera.SetPosition(vtkm::Vec3f_32(-2, 1.75, 1.75));
+  vtkm::cont::ColorTable colorTable("inferno");
+
+  // Background color:
+  vtkm::rendering::Color bg(0.2f, 0.2f, 0.2f, 1.0f);
+  vtkm::rendering::Scene scene;
+  int width = 512, height = 512;
+  CanvasRayTracer canvas(width, height);
+
+  for (int i = 0; i < dsPerRank; i++)
+  {
+    //Create a sequence of datasets along the X direction.
+    std::string fieldName = "tangle";
+    vtkm::source::Tangle tangle;
+    vtkm::Vec3f pt(rank * dsPerRank + i, 0, 0);
+    if (rank == 1)
+      std::cout << "PT= " << pt << std::endl;
+    tangle.SetPointDimensions({ 50, 50, 50 });
+    tangle.SetOrigin(pt);
+    vtkm::cont::DataSet ds = tangle.Execute();
+
+    vtkm::rendering::Actor actor(
+      ds.GetCellSet(), ds.GetCoordinateSystem(), ds.GetField(fieldName), colorTable);
+    scene.AddActor(actor);
+  }
+
+  vtkm::rendering::View3D view(scene, MapperRayTracer(), canvas, camera, bg);
+  view.Paint();
+
+  canvas.SaveAs("result.png");
+
+  /*
+  auto colors = &GetVTKMPointer(canvas.GetColorBuffer())[0][0];
+  auto depths = GetVTKMPointer(canvas.GetDepthBuffer());
+
+  vtkm::rendering::compositing::Compositor compositor;
+  compositor.AddImage(colors, depths, width, height);
+  auto res = compositor.Composite();
+  if (comm.rank() == 0)
+  {
+    res.Save("RESULT.png", { "" });
+  }
+  */
+}
+
+void TestVolumeRenderComposite(bool unstructured)
+{
+  using vtkm::rendering::CanvasRayTracer;
+  using vtkm::rendering::MapperConnectivity;
+  using vtkm::rendering::MapperVolume;
+
+  auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
+
+  int numBlocks = comm.size() * 1;
+  int rank = comm.rank();
+  int dsPerRank = 1;
+
+  vtkm::rendering::Camera camera;
+  camera.SetLookAt(vtkm::Vec3f_32(1.0, 0.5, 0.5));
+  camera.SetViewUp(vtkm::make_Vec(0.f, 1.f, 0.f));
+  camera.SetClippingRange(1.f, 10.f);
+  camera.SetFieldOfView(60.f);
+  camera.SetPosition(vtkm::Vec3f_32(-2, 1.75, 1.75));
+  vtkm::cont::ColorTable colorTable("inferno");
+
+  colorTable.AddPointAlpha(0.0, .01f);
+  colorTable.AddPointAlpha(1.0, .01f);
+
+  // Background color:
+  vtkm::rendering::Color bg(0.2f, 0.2f, 0.2f, 1.0f);
+  vtkm::rendering::Scene scene;
+  int width = 512, height = 512;
+  CanvasRayTracer canvas(width, height);
+
+  for (int i = 0; i < dsPerRank; i++)
+  {
+    //Create a sequence of datasets along the X direction.
+    std::string fieldName = "tangle";
+    vtkm::source::Tangle tangle;
+    vtkm::Vec3f pt(rank * dsPerRank + i, 0, 0);
+    if (rank == 1)
+      std::cout << "PT= " << pt << std::endl;
+    tangle.SetPointDimensions({ 50, 50, 50 });
+    tangle.SetOrigin(pt);
+    vtkm::cont::DataSet ds = tangle.Execute();
+
+    if (unstructured)
+    {
+      vtkm::filter::clean_grid::CleanGrid cleanGrid;
+      ds = cleanGrid.Execute(ds);
+    }
+
+    vtkm::rendering::Actor actor(
+      ds.GetCellSet(), ds.GetCoordinateSystem(), ds.GetField(fieldName), colorTable);
+    scene.AddActor(actor);
+  }
+
+
+
+  if (unstructured)
+  {
+    vtkm::rendering::View3D view(scene, MapperConnectivity(), canvas, camera, bg);
+    view.Paint();
+    canvas.SaveAs("result-unstructured.png");
+  }
+  else
+  {
+    vtkm::rendering::View3D view(scene, MapperVolume(), canvas, camera, bg);
+    view.Paint();
+    canvas.SaveAs("result-structured.png");
+  }
+}
+
+void RenderTests()
+{
+  //  TestImageComposite();
+  //TestRenderComposite();
+  TestVolumeRenderComposite(false);
+  TestVolumeRenderComposite(true);
+}
+
+} //namespace
+
+int UnitTestImageCompositing(int argc, char* argv[])
+{
+  return vtkm::cont::testing::Testing::Run(RenderTests, argc, argv);
+}

vtkm/rendering/testing/UnitTestMapperRayTracer.cxx

@@ -17,6 +17,9 @@
 #include <vtkm/rendering/View3D.h>
 #include <vtkm/rendering/testing/RenderTest.h>
 
+
+#include <vtkm/filter/field_transform/PointTransform.h>
+
 namespace
 {
 
@@ -29,8 +32,23 @@ void RenderTests()
   options.AllowAnyDevice = false;
   options.ColorTable = vtkm::cont::ColorTable::Preset::Inferno;
 
+
+  //  vtkm::rendering::testing::RenderTest(
+  //    maker.Make3DRegularDataSet0(), "pointvar", "rendering/raytracer/regular3D.png", options);
+
+  auto ds0 = maker.Make3DRegularDataSet0();
+  auto ds1 = maker.Make3DRegularDataSet0();
+
+  vtkm::filter::field_transform::PointTransform filter;
+  filter.SetTranslation({ -1, 2, 2 });
+  filter.SetOutputFieldName("coordinates");
+  filter.SetChangeCoordinateSystem(true);
+  auto res = filter.Execute(ds1);
+
   vtkm::rendering::testing::RenderTest(
-    maker.Make3DRegularDataSet0(), "pointvar", "rendering/raytracer/regular3D.png", options);
+    { { ds0, "pointvar" }, { res, "pointvar" } }, "rendering/raytracer/regular3D.png", options);
+
+
   vtkm::rendering::testing::RenderTest(maker.Make3DRectilinearDataSet0(),
                                        "pointvar",
                                        "rendering/raytracer/rectilinear3D.png",

vtkm/rendering/testing/t_vtkm_test_utils.h

@@ -0,0 +1,773 @@
+#ifndef t_test_utils_h
+#define t_test_utils_h
+
+#include <assert.h>
+#include <random>
+#include <vtkm/Matrix.h>
+#include <vtkm/VectorAnalysis.h>
+#include <vtkm/cont/DataSet.h>
+#include <vtkm/cont/DataSetBuilderExplicit.h>
+#include <vtkm/cont/DataSetBuilderRectilinear.h>
+#include <vtkm/cont/DataSetBuilderUniform.h>
+//#include <vtkm/cont/testing/Testing.h>
+
+#define BASE_SIZE 32
+typedef vtkm::cont::ArrayHandleUniformPointCoordinates UniformCoords;
+
+struct SpatialDivision
+{
+  int m_mins[3];
+  int m_maxs[3];
+
+  SpatialDivision()
+    : m_mins{ 0, 0, 0 }
+    , m_maxs{ 1, 1, 1 }
+  {
+  }
+
+  bool CanSplit(int dim) { return m_maxs[dim] - m_mins[dim] + 1 > 1; }
+
+  SpatialDivision Split(int dim)
+  {
+    SpatialDivision r_split;
+    r_split = *this;
+    assert(CanSplit(dim));
+    int size = m_maxs[dim] - m_mins[dim] + 1;
+    int left_offset = size / 2;
+
+    //shrink the left side
+    m_maxs[dim] = m_mins[dim] + left_offset - 1;
+    //shrink the right side
+    r_split.m_mins[dim] = m_maxs[dim] + 1;
+    return r_split;
+  }
+};
+
+SpatialDivision GetBlock(int block, int num_blocks, SpatialDivision total_size)
+{
+
+  std::vector<SpatialDivision> divs;
+  divs.push_back(total_size);
+  int avail = num_blocks - 1;
+  int current_dim = 0;
+  int missed_splits = 0;
+  const int num_dims = 3;
+  while (avail > 0)
+  {
+    const int current_size = divs.size();
+    int temp_avail = avail;
+    for (int i = 0; i < current_size; ++i)
+    {
+      if (avail == 0)
+        break;
+      if (!divs[i].CanSplit(current_dim))
+      {
+        continue;
+      }
+      divs.push_back(divs[i].Split(current_dim));
+      --avail;
+    }
+    if (temp_avail == avail)
+    {
+      // dims were too small to make any spit
+      missed_splits++;
+      if (missed_splits == 3)
+      {
+        // we tried all three dims and could
+        // not make a split.
+        for (int i = 0; i < avail; ++i)
+        {
+          SpatialDivision empty;
+          empty.m_maxs[0] = 0;
+          empty.m_maxs[1] = 0;
+          empty.m_maxs[2] = 0;
+          divs.push_back(empty);
+        }
+        if (block == 0)
+        {
+          std::cerr << "** Warning **: data set size is too small to"
+                    << " divide between " << num_blocks << " blocks. "
+                    << " Adding " << avail << " empty data sets\n";
+        }
+
+        avail = 0;
+      }
+    }
+    else
+    {
+      missed_splits = 0;
+    }
+
+    current_dim = (current_dim + 1) % num_dims;
+  }
+
+  return divs.at(block);
+}
+
+template <typename FieldType>
+vtkm::cont::Field CreateCellScalarField(int size, const char* fieldName)
+{
+  vtkm::cont::ArrayHandle<FieldType> data;
+  data.Allocate(size);
+
+  for (int i = 0; i < size; ++i)
+  {
+    FieldType val = i / vtkm::Float32(size);
+    data.WritePortal().Set(i, val);
+  }
+
+
+  vtkm::cont::Field field(fieldName, vtkm::cont::Field::Association::Cells, data);
+  return field;
+}
+
+vtkm::cont::Field CreateGhostScalarField(vtkm::Id3 dims)
+{
+  vtkm::Int32 size = dims[0] * dims[1] * dims[2];
+  vtkm::cont::ArrayHandle<vtkm::Int32> data;
+  data.Allocate(size);
+
+  for (int z = 0; z < dims[2]; ++z)
+    for (int y = 0; y < dims[1]; ++y)
+      for (int x = 0; x < dims[0]; ++x)
+      {
+        vtkm::UInt8 flag = 0;
+        if (x < 1 || x > dims[0] - 2)
+          flag = 1;
+        if (y < 1 || y > dims[1] - 2)
+          flag = 1;
+        if (z < 1 || z > dims[2] - 2)
+          flag = 1;
+        vtkm::Id index = z * dims[0] * dims[1] + y * dims[0] + x;
+        data.WritePortal().Set(index, flag);
+      }
+
+  vtkm::cont::Field field("ghosts", vtkm::cont::Field::Association::Cells, data);
+  return field;
+}
+
+template <typename FieldType>
+vtkm::cont::Field CreatePointScalarField(UniformCoords coords, const char* fieldName)
+
+{
+  const int size = coords.GetNumberOfValues();
+  vtkm::cont::ArrayHandle<FieldType> data;
+  data.Allocate(size);
+  auto portal = coords.ReadPortal();
+  for (int i = 0; i < size; ++i)
+  {
+    vtkm::Vec<FieldType, 3> point = portal.Get(i);
+
+    FieldType val = vtkm::Magnitude(point) + 1.f;
+    data.WritePortal().Set(i, val);
+  }
+
+  vtkm::cont::Field field(fieldName, vtkm::cont::Field::Association::Points, data);
+  return field;
+}
+
+template <typename FieldType>
+vtkm::cont::Field CreatePointVecField(int size, const char* fieldName)
+{
+  vtkm::cont::ArrayHandle<vtkm::Vec<FieldType, 3>> data;
+  data.Allocate(size);
+
+  for (int i = 0; i < size; ++i)
+  {
+    FieldType val = i / FieldType(size);
+
+    vtkm::Vec<FieldType, 3> vec(val, -val, val);
+
+    data.WritePortal().Set(i, vec);
+  }
+
+  vtkm::cont::Field field(fieldName, vtkm::cont::Field::Association::Points, data);
+  return field;
+}
+
+vtkm::cont::DataSet CreateTestData(int block, int num_blocks, int base_size)
+{
+  SpatialDivision mesh_size;
+
+  mesh_size.m_mins[0] = 0;
+  mesh_size.m_mins[1] = 0;
+  mesh_size.m_mins[2] = 0;
+
+  mesh_size.m_maxs[0] = num_blocks * base_size - 1;
+  mesh_size.m_maxs[1] = num_blocks * base_size - 1;
+  mesh_size.m_maxs[2] = num_blocks * base_size - 1;
+
+  SpatialDivision local_block = GetBlock(block, num_blocks, mesh_size);
+
+  vtkm::Vec<vtkm::Float32, 3> origin;
+  origin[0] = local_block.m_mins[0];
+  origin[1] = local_block.m_mins[1];
+  origin[2] = local_block.m_mins[2];
+
+  vtkm::Vec<vtkm::Float32, 3> spacing(1.f, 1.f, 1.f);
+
+  vtkm::Id3 point_dims;
+  point_dims[0] = local_block.m_maxs[0] - local_block.m_mins[0] + 2;
+  point_dims[1] = local_block.m_maxs[1] - local_block.m_mins[1] + 2;
+  point_dims[2] = local_block.m_maxs[2] - local_block.m_mins[2] + 2;
+
+
+  vtkm::Id3 cell_dims;
+  cell_dims[0] = point_dims[0] - 1;
+  cell_dims[1] = point_dims[1] - 1;
+  cell_dims[2] = point_dims[2] - 1;
+
+  vtkm::cont::DataSet data_set;
+
+  UniformCoords point_handle(point_dims, origin, spacing);
+
+  vtkm::cont::CoordinateSystem coords("coords", point_handle);
+  data_set.AddCoordinateSystem(coords);
+
+  vtkm::cont::CellSetStructured<3> cell_set;
+  cell_set.SetPointDimensions(point_dims);
+  data_set.SetCellSet(cell_set);
+
+  int num_points = point_dims[0] * point_dims[1] * point_dims[2];
+  int num_cells = cell_dims[0] * cell_dims[1] * cell_dims[2];
+
+  data_set.AddField(CreatePointScalarField<vtkm::Float32>(point_handle, "point_data_Float32"));
+  data_set.AddField(CreatePointVecField<vtkm::Float32>(num_points, "vector_data_Float32"));
+  data_set.AddField(CreateCellScalarField<vtkm::Float32>(num_cells, "cell_data_Float32"));
+  data_set.AddField(CreatePointScalarField<vtkm::Float64>(point_handle, "point_data_Float64"));
+  data_set.AddField(CreatePointVecField<vtkm::Float64>(num_points, "vector_data_Float64"));
+  data_set.AddField(CreateCellScalarField<vtkm::Float64>(num_cells, "cell_data_Float64"));
+  data_set.AddField(CreateGhostScalarField(cell_dims));
+  return data_set;
+}
+
+vtkm::cont::DataSet CreateTestDataRectilinear(int block, int num_blocks, int base_size)
+{
+  SpatialDivision mesh_size;
+
+  mesh_size.m_mins[0] = 0;
+  mesh_size.m_mins[1] = 0;
+  mesh_size.m_mins[2] = 0;
+
+  mesh_size.m_maxs[0] = num_blocks * base_size - 1;
+  mesh_size.m_maxs[1] = num_blocks * base_size - 1;
+  mesh_size.m_maxs[2] = num_blocks * base_size - 1;
+
+  SpatialDivision local_block = GetBlock(block, num_blocks, mesh_size);
+
+  vtkm::Vec<vtkm::Float32, 3> origin;
+  origin[0] = local_block.m_mins[0];
+  origin[1] = local_block.m_mins[1];
+  origin[2] = local_block.m_mins[2];
+
+  vtkm::Vec<vtkm::Float32, 3> spacing(1.f, 1.f, 1.f);
+
+  vtkm::Id3 point_dims;
+  point_dims[0] = local_block.m_maxs[0] - local_block.m_mins[0] + 2;
+  point_dims[1] = local_block.m_maxs[1] - local_block.m_mins[1] + 2;
+  point_dims[2] = local_block.m_maxs[2] - local_block.m_mins[2] + 2;
+
+
+  vtkm::Id3 cell_dims;
+  cell_dims[0] = point_dims[0] - 1;
+  cell_dims[1] = point_dims[1] - 1;
+  cell_dims[2] = point_dims[2] - 1;
+
+  std::vector<vtkm::Float64> xvals, yvals, zvals;
+  xvals.resize((size_t)point_dims[0]);
+  xvals[0] = static_cast<vtkm::Float64>(local_block.m_mins[0]);
+  for (size_t i = 1; i < (size_t)point_dims[0]; i++)
+    xvals[i] = xvals[i - 1] + spacing[0];
+
+  yvals.resize((size_t)point_dims[1]);
+  yvals[0] = static_cast<vtkm::Float64>(local_block.m_mins[1]);
+  for (size_t i = 1; i < (size_t)point_dims[1]; i++)
+    yvals[i] = yvals[i - 1] + spacing[1];
+
+  zvals.resize((size_t)point_dims[2]);
+  zvals[0] = static_cast<vtkm::Float64>(local_block.m_mins[2]);
+  for (size_t i = 1; i < (size_t)point_dims[2]; i++)
+    zvals[i] = zvals[i - 1] + spacing[2];
+
+  vtkm::cont::DataSetBuilderRectilinear dataSetBuilder;
+  vtkm::cont::DataSet data_set = dataSetBuilder.Create(xvals, yvals, zvals);
+
+  int num_points = point_dims[0] * point_dims[1] * point_dims[2];
+
+  data_set.AddField(CreatePointVecField<vtkm::Float32>(num_points, "vector_data_Float32"));
+  data_set.AddField(CreatePointVecField<vtkm::Float64>(num_points, "vector_data_Float64"));
+
+  return data_set;
+}
+
+vtkm::cont::DataSet CreateTestDataPoints(int num_points)
+{
+  std::vector<double> x_vals;
+  std::vector<double> y_vals;
+  std::vector<double> z_vals;
+  std::vector<vtkm::UInt8> shapes;
+  std::vector<vtkm::IdComponent> num_indices;
+  std::vector<vtkm::Id> conn;
+  std::vector<double> field;
+
+  x_vals.resize(num_points);
+  y_vals.resize(num_points);
+  z_vals.resize(num_points);
+  shapes.resize(num_points);
+  conn.resize(num_points);
+  num_indices.resize(num_points);
+  field.resize(num_points);
+
+  std::linear_congruential_engine<std::uint_fast32_t, 48271, 0, 2147483647> rgen{ 0 };
+  std::uniform_real_distribution<double> dist{ -10., 10. };
+
+  for (int i = 0; i < num_points; ++i)
+  {
+    x_vals[i] = dist(rgen);
+    y_vals[i] = dist(rgen);
+    z_vals[i] = dist(rgen);
+    field[i] = dist(rgen);
+    shapes[i] = vtkm::CELL_SHAPE_VERTEX;
+    num_indices[i] = 1;
+    conn[i] = i;
+  }
+  vtkm::cont::DataSetBuilderExplicit dataSetBuilder;
+  vtkm::cont::DataSet data_set =
+    dataSetBuilder.Create(x_vals, y_vals, z_vals, shapes, num_indices, conn);
+  vtkm::cont::Field vfield = vtkm::cont::make_Field(
+    "point_data_Float64", vtkm::cont::Field::Association::Points, field, vtkm::CopyFlag::On);
+  data_set.AddField(vfield);
+  return data_set;
+}
+
+//-----------------------------------------------------------------------------
+//Create VTK-m Data Sets
+//-----------------------------------------------------------------------------
+
+//Make a 2Duniform dataset.
+inline vtkm::cont::DataSet Make2DUniformDataSet0()
+{
+  vtkm::cont::DataSetBuilderUniform dsb;
+  constexpr vtkm::Id2 dimensions(3, 2);
+  vtkm::cont::DataSet dataSet = dsb.Create(dimensions);
+
+  constexpr vtkm::Id nVerts = 6;
+  constexpr vtkm::Float32 var[nVerts] = { 10.1f, 20.1f, 30.1f, 40.1f, 50.1f, 60.1f };
+
+  dataSet.AddPointField("pointvar", var, nVerts);
+
+  constexpr vtkm::Float32 cellvar[2] = { 100.1f, 200.1f };
+  dataSet.AddCellField("cellvar", cellvar, 2);
+
+  return dataSet;
+}
+
+//Make a 2D rectilinear dataset.
+inline vtkm::cont::DataSet Make2DRectilinearDataSet0()
+{
+  vtkm::cont::DataSetBuilderRectilinear dsb;
+  std::vector<vtkm::Float32> X(3), Y(2);
+
+  X[0] = 0.0f;
+  X[1] = 1.0f;
+  X[2] = 2.0f;
+  Y[0] = 0.0f;
+  Y[1] = 1.0f;
+
+  vtkm::cont::DataSet dataSet = dsb.Create(X, Y);
+
+  const vtkm::Id nVerts = 6;
+  vtkm::Float32 var[nVerts];
+  for (int i = 0; i < nVerts; i++)
+    var[i] = (vtkm::Float32)i;
+  dataSet.AddPointField("pointvar", var, nVerts);
+
+  const vtkm::Id nCells = 2;
+  vtkm::Float32 cellvar[nCells];
+  for (int i = 0; i < nCells; i++)
+    cellvar[i] = (vtkm::Float32)i;
+  dataSet.AddCellField("cellvar", cellvar, nCells);
+
+  return dataSet;
+}
+
+inline vtkm::cont::DataSet Make3DExplicitDataSet5()
+{
+  vtkm::cont::DataSet dataSet;
+
+  const int nVerts = 11;
+  using CoordType = vtkm::Vec3f_32;
+  CoordType coordinates[nVerts] = {
+    CoordType(0, 0, 0),     //0
+    CoordType(1, 0, 0),     //1
+    CoordType(1, 0, 1),     //2
+    CoordType(0, 0, 1),     //3
+    CoordType(0, 1, 0),     //4
+    CoordType(1, 1, 0),     //5
+    CoordType(1, 1, 1),     //6
+    CoordType(0, 1, 1),     //7
+    CoordType(2, 0.5, 0.5), //8
+    CoordType(0, 2, 0),     //9
+    CoordType(1, 2, 0)      //10
+  };
+  vtkm::Float32 vars[nVerts] = { 10.1f, 20.1f, 30.2f, 40.2f, 50.3f, 60.2f,
+                                 70.2f, 80.3f, 90.f,  10.f,  11.f };
+
+  dataSet.AddCoordinateSystem(
+    vtkm::cont::make_CoordinateSystem("coordinates", coordinates, nVerts, vtkm::CopyFlag::On));
+
+  //Set point scalar
+  dataSet.AddField(make_Field(
+    "pointvar", vtkm::cont::Field::Association::Points, vars, nVerts, vtkm::CopyFlag::On));
+
+  //Set cell scalar
+  const int nCells = 4;
+  vtkm::Float32 cellvar[nCells] = { 100.1f, 110.f, 120.2f, 130.5f };
+  dataSet.AddField(make_Field(
+    "cellvar", vtkm::cont::Field::Association::Cells, cellvar, nCells, vtkm::CopyFlag::On));
+
+  vtkm::cont::CellSetExplicit<> cellSet;
+  vtkm::Vec<vtkm::Id, 8> ids;
+
+  cellSet.PrepareToAddCells(nCells, 23);
+
+  ids[0] = 0;
+  ids[1] = 1;
+  ids[2] = 5;
+  ids[3] = 4;
+  ids[4] = 3;
+  ids[5] = 2;
+  ids[6] = 6;
+  ids[7] = 7;
+  cellSet.AddCell(vtkm::CELL_SHAPE_HEXAHEDRON, 8, ids);
+
+  ids[0] = 1;
+  ids[1] = 5;
+  ids[2] = 6;
+  ids[3] = 2;
+  ids[4] = 8;
+  cellSet.AddCell(vtkm::CELL_SHAPE_PYRAMID, 5, ids);
+
+  ids[0] = 5;
+  ids[1] = 8;
+  ids[2] = 10;
+  ids[3] = 6;
+  cellSet.AddCell(vtkm::CELL_SHAPE_TETRA, 4, ids);
+
+  ids[0] = 4;
+  ids[1] = 7;
+  ids[2] = 9;
+  ids[3] = 5;
+  ids[4] = 6;
+  ids[5] = 10;
+  cellSet.AddCell(vtkm::CELL_SHAPE_WEDGE, 6, ids);
+
+  cellSet.CompleteAddingCells(nVerts);
+
+  //todo this need to be a reference/shared_ptr style class
+  dataSet.SetCellSet(cellSet);
+
+  return dataSet;
+}
+
+inline vtkm::cont::DataSet Make3DUniformDataSet0()
+{
+  vtkm::cont::DataSetBuilderUniform dsb;
+  constexpr vtkm::Id3 dimensions(3, 2, 3);
+  vtkm::cont::DataSet dataSet = dsb.Create(dimensions);
+
+  constexpr int nVerts = 18;
+  constexpr vtkm::Float32 vars[nVerts] = { 10.1f,  20.1f,  30.1f,  40.1f,  50.2f,  60.2f,
+                                           70.2f,  80.2f,  90.3f,  100.3f, 110.3f, 120.3f,
+                                           130.4f, 140.4f, 150.4f, 160.4f, 170.5f, 180.5f };
+
+  //Set point and cell scalar
+  dataSet.AddPointField("pointvar", vars, nVerts);
+
+  constexpr vtkm::Float32 cellvar[4] = { 100.1f, 100.2f, 100.3f, 100.4f };
+  dataSet.AddCellField("cellvar", cellvar, 4);
+
+  return dataSet;
+}
+
+inline vtkm::cont::DataSet Make3DExplicitDataSet2()
+{
+  vtkm::cont::DataSet dataSet;
+
+  const int nVerts = 8;
+  using CoordType = vtkm::Vec3f_32;
+  CoordType coordinates[nVerts] = {
+    CoordType(0, 0, 0), // 0
+    CoordType(1, 0, 0), // 1
+    CoordType(1, 0, 1), // 2
+    CoordType(0, 0, 1), // 3
+    CoordType(0, 1, 0), // 4
+    CoordType(1, 1, 0), // 5
+    CoordType(1, 1, 1), // 6
+    CoordType(0, 1, 1)  // 7
+  };
+  vtkm::Float32 vars[nVerts] = { 10.1f, 20.1f, 30.2f, 40.2f, 50.3f, 60.2f, 70.2f, 80.3f };
+
+  dataSet.AddCoordinateSystem(
+    vtkm::cont::make_CoordinateSystem("coordinates", coordinates, nVerts, vtkm::CopyFlag::On));
+
+  //Set point scalar
+  dataSet.AddField(make_Field(
+    "pointvar", vtkm::cont::Field::Association::Points, vars, nVerts, vtkm::CopyFlag::On));
+
+  //Set cell scalar
+  vtkm::Float32 cellvar[2] = { 100.1f };
+  dataSet.AddField(
+    make_Field("cellvar", vtkm::cont::Field::Association::Cells, cellvar, 1, vtkm::CopyFlag::On));
+
+  vtkm::cont::CellSetExplicit<> cellSet;
+  vtkm::Vec<vtkm::Id, 8> ids;
+  ids[0] = 0;
+  ids[1] = 1;
+  ids[2] = 2;
+  ids[3] = 3;
+  ids[4] = 4;
+  ids[5] = 5;
+  ids[6] = 6;
+  ids[7] = 7;
+
+  cellSet.PrepareToAddCells(1, 8);
+  cellSet.AddCell(vtkm::CELL_SHAPE_HEXAHEDRON, 8, ids);
+  cellSet.CompleteAddingCells(nVerts);
+
+  //todo this need to be a reference/shared_ptr style class
+  dataSet.SetCellSet(cellSet);
+
+  return dataSet;
+}
+
+#if 0
+
+namespace detail
+{
+
+template <typename T>
+struct TestValueImpl;
+} //namespace detail
+
+// Many tests involve getting and setting values in some index-based structure
+// (like an array). These tests also often involve trying many types. The
+// overloaded TestValue function returns some unique value for an index for a
+// given type. Different types might give different values.
+//
+template <typename T>
+static inline T TestValue(vtkm::Id index, T)
+{
+  return detail::TestValueImpl<T>()(index);
+}
+
+namespace detail
+{
+
+template <typename T>
+struct TestValueImpl
+{
+  T DoIt(vtkm::Id index, vtkm::TypeTraitsIntegerTag) const
+  {
+    constexpr bool larger_than_2bytes = sizeof(T) > 2;
+    if (larger_than_2bytes)
+    {
+       return T(index * 100);
+    }
+    else
+    {
+      return T(index + 100);
+    }
+  }
+
+  T DoIt(vtkm::Id index, vtkm::TypeTraitsRealTag) const
+  {
+    return T(0.01f * static_cast<float>(index) + 1.001f);
+  }
+
+  T operator()(vtkm::Id index) const
+  {
+    return this->DoIt(index, typename vtkm::TypeTraits<T>::NumericTag());
+  }
+};
+
+template <typename T, vtkm::IdComponent N>
+struct TestValueImpl<vtkm::Vec<T, N>>
+{
+  vtkm::Vec<T, N> operator()(vtkm::Id index) const
+  {
+    vtkm::Vec<T, N> value;
+    for (vtkm::IdComponent i = 0; i < N; i++)
+    {
+      value[i] = TestValue(index * N + i, T());
+    }
+    return value;
+  }
+};
+
+template <typename U, typename V>
+struct TestValueImpl<vtkm::Pair<U, V>>
+{
+  vtkm::Pair<U, V> operator()(vtkm::Id index) const
+  {
+    return vtkm::Pair<U, V>(TestValue(2 * index, U()), TestValue(2 * index + 1, V()));
+  }
+};
+
+template <typename T, vtkm::IdComponent NumRow, vtkm::IdComponent NumCol>
+struct TestValueImpl<vtkm::Matrix<T, NumRow, NumCol>>
+{
+  vtkm::Matrix<T, NumRow, NumCol> operator()(vtkm::Id index) const
+  {
+    vtkm::Matrix<T, NumRow, NumCol> value;
+    vtkm::Id runningIndex = index * NumRow * NumCol;
+    for (vtkm::IdComponent row = 0; row < NumRow; ++row)
+    {
+      for (vtkm::IdComponent col = 0; col < NumCol; ++col)
+      {
+        value(row, col) = TestValue(runningIndex, T());
+        ++runningIndex;
+      }
+    }
+    return value;
+  }
+};
+
+template <>
+struct TestValueImpl<std::string>
+{
+  std::string operator()(vtkm::Id index) const
+  {
+    std::stringstream stream;
+    stream << index;
+    return stream.str();
+  }
+};
+
+} //namespace detail
+
+// Verifies that the contents of the given array portal match the values
+// returned by vtkm::testing::TestValue.
+template <typename PortalType>
+static inline void CheckPortal(const PortalType& portal)
+{
+  using ValueType = typename PortalType::ValueType;
+  for (vtkm::Id index = 0; index < portal.GetNumberOfValues(); index++)
+  {
+    ValueType expectedValue = TestValue(index, ValueType());
+    ValueType foundValue = portal.Get(index);
+    if (!test_equal(expectedValue, foundValue))
+    {
+      ASCENT_ERROR("Got unexpected value in array. Expected: " << expectedValue
+              << ", Found: " << foundValue << "\n");
+    }
+  }
+}
+
+/// Sets all the values in a given array portal to be the values returned
+/// by vtkm::testing::TestValue. The ArrayPortal must be allocated first.
+template <typename PortalType>
+static inline void SetPortal(const PortalType& portal)
+{
+  using ValueType = typename PortalType::ValueType;
+
+  for (vtkm::Id index = 0; index < portal.GetNumberOfValues(); index++)
+  {
+    portal.Set(index, TestValue(index, ValueType()));
+  }
+}
+
+#endif
+
+inline vtkm::cont::DataSet Make3DExplicitDataSetCowNose()
+{
+  // prepare data array
+  const int nVerts = 17;
+  using CoordType = vtkm::Vec3f_64;
+  CoordType coordinates[nVerts] = {
+    CoordType(0.0480879, 0.151874, 0.107334),     CoordType(0.0293568, 0.245532, 0.125337),
+    CoordType(0.0224398, 0.246495, 0.1351),       CoordType(0.0180085, 0.20436, 0.145316),
+    CoordType(0.0307091, 0.152142, 0.0539249),    CoordType(0.0270341, 0.242992, 0.107567),
+    CoordType(0.000684071, 0.00272505, 0.175648), CoordType(0.00946217, 0.077227, 0.187097),
+    CoordType(-0.000168991, 0.0692243, 0.200755), CoordType(-0.000129414, 0.00247137, 0.176561),
+    CoordType(0.0174172, 0.137124, 0.124553),     CoordType(0.00325994, 0.0797155, 0.184912),
+    CoordType(0.00191765, 0.00589327, 0.16608),   CoordType(0.0174716, 0.0501928, 0.0930275),
+    CoordType(0.0242103, 0.250062, 0.126256),     CoordType(0.0108188, 0.152774, 0.167914),
+    CoordType(5.41687e-05, 0.00137834, 0.175119)
+  };
+  const int connectivitySize = 57;
+  vtkm::Id pointId[connectivitySize] = { 0, 1, 3,  2, 3,  1, 4,  5,  0,  1, 0,  5,  7,  8,  6,
+                                         9, 6, 8,  0, 10, 7, 11, 7,  10, 0, 6,  13, 12, 13, 6,
+                                         1, 5, 14, 1, 14, 2, 0,  3,  15, 0, 13, 4,  6,  16, 12,
+                                         6, 9, 16, 7, 11, 8, 0,  15, 10, 7, 6,  0 };
+
+  // create DataSet
+  vtkm::cont::DataSet dataSet;
+  dataSet.AddCoordinateSystem(
+    vtkm::cont::make_CoordinateSystem("coordinates", coordinates, nVerts, vtkm::CopyFlag::On));
+
+  vtkm::cont::ArrayHandle<vtkm::Id> connectivity;
+  connectivity.Allocate(connectivitySize);
+
+  for (vtkm::Id i = 0; i < connectivitySize; ++i)
+  {
+    connectivity.WritePortal().Set(i, pointId[i]);
+  }
+  vtkm::cont::CellSetSingleType<> cellSet;
+  cellSet.Fill(nVerts, vtkm::CELL_SHAPE_TRIANGLE, 3, connectivity);
+  dataSet.SetCellSet(cellSet);
+
+  std::vector<vtkm::Float32> pointvar(nVerts);
+  std::iota(pointvar.begin(), pointvar.end(), 15.f);
+  std::vector<vtkm::Float32> cellvar(connectivitySize / 3);
+  std::iota(cellvar.begin(), cellvar.end(), 132.f);
+
+  vtkm::cont::ArrayHandle<vtkm::Vec3f> pointvec;
+  pointvec.Allocate(nVerts);
+  SetPortal(pointvec.WritePortal());
+
+  vtkm::cont::ArrayHandle<vtkm::Vec3f> cellvec;
+  cellvec.Allocate(connectivitySize / 3);
+  SetPortal(cellvec.WritePortal());
+
+  dataSet.AddPointField("pointvar", pointvar);
+  dataSet.AddCellField("cellvar", cellvar);
+  dataSet.AddPointField("point_vectors", pointvec);
+  dataSet.AddCellField("cell_vectors", cellvec);
+
+  return dataSet;
+}
+
+inline vtkm::cont::DataSet Make3DRectilinearDataSet0()
+{
+  vtkm::cont::DataSetBuilderRectilinear dsb;
+  std::vector<vtkm::Float32> X(3), Y(2), Z(3);
+
+  X[0] = 0.0f;
+  X[1] = 1.0f;
+  X[2] = 2.0f;
+  Y[0] = 0.0f;
+  Y[1] = 1.0f;
+  Z[0] = 0.0f;
+  Z[1] = 1.0f;
+  Z[2] = 2.0f;
+
+  vtkm::cont::DataSet dataSet = dsb.Create(X, Y, Z);
+
+  const vtkm::Id nVerts = 18;
+  vtkm::Float32 var[nVerts];
+  for (int i = 0; i < nVerts; i++)
+    var[i] = (vtkm::Float32)i;
+  dataSet.AddPointField("pointvar", var, nVerts);
+
+  const vtkm::Id nCells = 4;
+  vtkm::Float32 cellvar[nCells];
+  for (int i = 0; i < nCells; i++)
+    cellvar[i] = (vtkm::Float32)i;
+  dataSet.AddCellField("cellvar", cellvar, nCells);
+
+  return dataSet;
+}
+
+
+#endif

vtkm/rendering/testlib/RenderTest.cxx

@@ -173,8 +173,11 @@ void DoRenderTest(vtkm::rendering::Canvas& canvas,
   vtkm::Range fieldRange;
   for (std::size_t dataFieldId = 0; dataFieldId < numFields; ++dataFieldId)
   {
+    std::cout << "   AddActor: " << dataFieldId << std::endl;
+
     vtkm::cont::DataSet dataSet = dataSetsFields[dataFieldId].first;
     std::string fieldName = dataSetsFields[dataFieldId].second;
+    dataSet.PrintSummary(std::cout);
     if (options.Colors.empty())
     {
       scene.AddActor(vtkm::rendering::Actor(dataSet.GetCellSet(),
@@ -280,6 +283,7 @@ void RenderTest(const vtkm::cont::DataSet& dataSet,
                 const std::string& outputFile,
                 const RenderTestOptions& options)
 {
+  std::cout << "RenderTest_1!!!" << std::endl;
   RenderTest({ { dataSet, fieldName } }, outputFile, options);
 }
 
@@ -287,6 +291,7 @@ void RenderTest(const DataSetFieldVector& dataSetsFields,
                 const std::string& outputFile,
                 const RenderTestOptions& options)
 {
+  std::cout << "RenderTest_2!!!" << std::endl;
   std::unique_ptr<vtkm::cont::ScopedRuntimeDeviceTracker> deviceScope;
   if (options.AllowAnyDevice)
   {

vtkm/source/PerlinNoise.cxx

@@ -212,8 +212,8 @@ vtkm::cont::DataSet PerlinNoise::DoExecute() const
 
   vtkm::cont::DataSet dataSet;
   const vtkm::Vec3f cellDims = this->GetCellDimensions();
-  const vtkm::Vec3f spacing(1.0f / cellDims[0], 1.0f / cellDims[1], 1.0f / cellDims[2]);
-
+  const vtkm::Vec3f extentDelta = this->MaxExtent - this->Origin;
+  const vtkm::Vec3f spacing = extentDelta / cellDims;
 
   vtkm::cont::CellSetStructured<3> cellSet;
   cellSet.SetPointDimensions(this->PointDimensions);

vtkm/source/PerlinNoise.h

@@ -55,6 +55,9 @@ public:
   VTKM_CONT vtkm::Vec3f GetOrigin() const { return this->Origin; }
   VTKM_CONT void SetOrigin(const vtkm::Vec3f& origin) { this->Origin = origin; }
 
+  VTKM_CONT vtkm::Vec3f GetMaxExtent() const { return this->MaxExtent; }
+  VTKM_CONT void SetMaxExtent(const vtkm::Vec3f& maxExtent) { this->MaxExtent = maxExtent; }
+
   /// \brief The seed used for the pseudorandom number generation of the noise.
   ///
   /// If the seed is not set, then a new, unique seed is picked each time `Execute` is run.
@@ -70,6 +73,7 @@ private:
 
   vtkm::Id3 PointDimensions = { 16, 16, 16 };
   vtkm::Vec3f Origin = { 0, 0, 0 };
+  vtkm::Vec3f MaxExtent = { 1, 1, 1 };
   vtkm::IdComponent Seed = 0;
   bool SeedSet = false;
 };

vtkm/source/Tangle.cxx

@@ -78,13 +78,12 @@ vtkm::cont::DataSet Tangle::DoExecute() const
   vtkm::cont::ArrayHandle<vtkm::Float32> pointFieldArray;
   this->Invoke(tangle::TangleField{ cellDims, mins, maxs }, cellSet, pointFieldArray);
 
-  const vtkm::Vec3f origin(0.0f, 0.0f, 0.0f);
   const vtkm::Vec3f spacing(1.0f / static_cast<vtkm::FloatDefault>(cellDims[0]),
                             1.0f / static_cast<vtkm::FloatDefault>(cellDims[1]),
                             1.0f / static_cast<vtkm::FloatDefault>(cellDims[2]));
 
   vtkm::cont::ArrayHandleUniformPointCoordinates coordinates(
-    this->PointDimensions, origin, spacing);
+    this->PointDimensions, this->Origin, spacing);
   dataSet.AddCoordinateSystem(vtkm::cont::CoordinateSystem("coordinates", coordinates));
   dataSet.AddField(vtkm::cont::make_FieldPoint("tangle", pointFieldArray));
 

vtkm/source/Tangle.h

@@ -44,6 +44,9 @@ public:
   {
   }
 
+  VTKM_CONT vtkm::Vec3f GetOrigin() const { return this->Origin; }
+  VTKM_CONT void SetOrigin(vtkm::Vec3f& pt) { this->Origin = pt; }
+
   VTKM_CONT vtkm::Id3 GetPointDimensions() const { return this->PointDimensions; }
   VTKM_CONT void SetPointDimensions(vtkm::Id3 dims) { this->PointDimensions = dims; }
 
@@ -54,6 +57,7 @@ private:
   vtkm::cont::DataSet DoExecute() const override;
 
   vtkm::Id3 PointDimensions = { 16, 16, 16 };
+  vtkm::Vec3f Origin = { 0, 0, 0 };
 };
 } //namespace source
 } //namespace vtkm

vtkm/thirdparty/diy/CMakeLists.txt

@@ -120,11 +120,18 @@ vtkm_install_targets(TARGETS vtkm_diy)
 if (NOT VTKm_INSTALL_ONLY_LIBRARIES)
   install(FILES
     ${VTKm_BINARY_INCLUDE_DIR}/${kit_dir}/Configure.h
+    ${CMAKE_CURRENT_SOURCE_DIR}/assigner.h
+    ${CMAKE_CURRENT_SOURCE_DIR}/decomposition.h
     ${CMAKE_CURRENT_SOURCE_DIR}/diy.h
     ${CMAKE_CURRENT_SOURCE_DIR}/environment.h
+    ${CMAKE_CURRENT_SOURCE_DIR}/master.h
     ${CMAKE_CURRENT_SOURCE_DIR}/mpi-cast.h
+    ${CMAKE_CURRENT_SOURCE_DIR}/point.h
     ${CMAKE_CURRENT_SOURCE_DIR}/post-include.h
     ${CMAKE_CURRENT_SOURCE_DIR}/pre-include.h
+    ${CMAKE_CURRENT_SOURCE_DIR}/reduce.h
+    ${CMAKE_CURRENT_SOURCE_DIR}/reduce-operations.h
     ${CMAKE_CURRENT_SOURCE_DIR}/serialization.h
+    ${CMAKE_CURRENT_SOURCE_DIR}/swap.h
     DESTINATION ${VTKm_INSTALL_INCLUDE_DIR}/${kit_dir}/)
 endif()

vtkm/thirdparty/diy/assigner.h

@@ -0,0 +1,19 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+#ifndef vtk_m_thirdparty_diy_assigner_h
+#define vtk_m_thirdparty_diy_assigner_h
+
+#include "pre-include.h"
+// clang-format off
+#include VTKM_DIY_INCLUDE(assigner.hpp)
+// clang-format on
+#include "post-include.h"
+
+#endif // vtk_m_thirdparty_diy_assigner_h

vtkm/thirdparty/diy/decomposition.h

@@ -0,0 +1,19 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+#ifndef vtk_m_thirdparty_diy_decomposition_h
+#define vtk_m_thirdparty_diy_decomposition_h
+
+#include "pre-include.h"
+// clang-format off
+#include VTKM_DIY_INCLUDE(decomposition.hpp)
+// clang-format on
+#include "post-include.h"
+
+#endif // vtk_m_thirdparty_diy_decomposition_h

vtkm/thirdparty/diy/master.h

@@ -0,0 +1,19 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+#ifndef vtk_m_thirdparty_diy_master_h
+#define vtk_m_thirdparty_diy_master_h
+
+#include "pre-include.h"
+// clang-format off
+#include VTKM_DIY_INCLUDE(master.hpp)
+// clang-format on
+#include "post-include.h"
+
+#endif // vtk_m_thirdparty_diy_master_h

vtkm/thirdparty/diy/point.h

@@ -0,0 +1,19 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+#ifndef vtk_m_thirdparty_diy_point_h
+#define vtk_m_thirdparty_diy_point_h
+
+#include "pre-include.h"
+// clang-format off
+#include VTKM_DIY_INCLUDE(point.hpp)
+// clang-format on
+#include "post-include.h"
+
+#endif // vtk_m_thirdparty_diy_point_h

vtkm/thirdparty/diy/reduce-operations.h

@@ -0,0 +1,20 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_thirdparty_diy_reduce_operations_h
+#define vtk_m_thirdparty_diy_reduce_operations_h
+
+#include "pre-include.h"
+// clang-format off
+#include VTKM_DIY_INCLUDE(reduce-operations.hpp)
+// clang-format on
+#include "post-include.h"
+
+#endif // vtk_m_thirdparty_diy_reduce_operations_h

vtkm/thirdparty/diy/reduce.h

@@ -0,0 +1,20 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+
+#ifndef vtk_m_thirdparty_diy_reduce_h
+#define vtk_m_thirdparty_diy_reduce_h
+
+#include "pre-include.h"
+// clang-format off
+#include VTKM_DIY_INCLUDE(reduce.hpp)
+// clang-format on
+#include "post-include.h"
+
+#endif // vtk_m_thirdparty_diy_reduce_reduce_h

vtkm/thirdparty/diy/swap.h

@@ -0,0 +1,19 @@
+//============================================================================
+//  Copyright (c) Kitware, Inc.
+//  All rights reserved.
+//  See LICENSE.txt for details.
+//
+//  This software is distributed WITHOUT ANY WARRANTY; without even
+//  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+//  PURPOSE.  See the above copyright notice for more information.
+//============================================================================
+#ifndef vtk_m_thirdparty_diy_swap_h
+#define vtk_m_thirdparty_diy_swap_h
+
+#include "pre-include.h"
+// clang-format off
+#include VTKM_DIY_INCLUDE(partners/swap.hpp)
+// clang-format on
+#include "post-include.h"
+
+#endif // vtk_m_thirdparty_diy_swap_h