Commit a632f77e authored by Alvaro Sanchez's avatar Alvaro Sanchez

Added support for multiple volume inputs in GPUVolumeRayCastMapper

When multiple inputs are rendered through the same mapper, it is possible to
composite overlapping areas correctly. Inputs are connected directly to
the mapper and their parameters (transfer functions, transformations, etc.) are
specified through standard vtkVolume instances. These vtkVolume instances are
to be registered in a special vtkProp3D, vtkMultiVolume.  vtkMultiVolume
represents a 3D space containing a set of vtkVolume instances which could be
but are not necessarily overlapping, this instance computes the bounding-box
containing all of its registered vtkVolumes.

vtkGPUVolumeRayCastMapper prepares and keeps track of the inputs as they are
allowed to be disconnected at any time. vtkOpenGLGPUVRCMapper's internals
particular to an input (vtkVolumeTexture, OpenGL helper instances, etc.) have
been promoted to a separate class vtkVolumeInputHelper. Other internals (such as
ComputeBounds, ComputePointToCellMatrix and related variables) have been moved
into vtkVolumeTexture as they are inherent to the texture definition (or
some of its blocks during streaming).

A separate code path is used when rendering multiple-inputs in order to
facilitate the co-existance of these two modes (single/multiple), due to current
feature incompatibilities with multiple inputs (e.g. texture-streaming, cropping,
etc.). Splitting DoGPURender into various smaller routines improves code
reusability between these two paths.

A limited set of the mapper features are currently supported for multiple
inputs:

  * Transfer functions: these are defined separately for per input.
    + 1D color
    + 1D scalar opacity
    + 1D gradient magnitude opacity
    + 2D scalar-gradient magnitude

  * Point and cell data
    + With the limitation that all of the inputs are assumed to share the same
    name/id.

  * Composite blending (front-to-back)

Test cases:

  * TestGPURayCastMultiVolumeOverlapping
  * TestGPURayCastMultiVolumeTransfer2D
  * TestGPURayCastMultiVolumeAddRemove
  * TestGPURayCastMultiVolumeCellData
parent feb05ce0
......@@ -301,6 +301,87 @@ public:
};
//@}
// .NAME vtkVector4 - templated base type for storage of 4D vectors.
//
template<typename T>
class vtkVector4 : public vtkVector<T, 4>
{
public:
vtkVector4()
{
}
explicit vtkVector4(const T& scalar) : vtkVector<T, 4>(scalar)
{
}
explicit vtkVector4(const T* init) : vtkVector<T, 4>(init)
{
}
vtkVector4(const T& x, const T& y, const T& z, const T& w)
{
this->Data[0] = x;
this->Data[1] = y;
this->Data[2] = z;
this->Data[3] = w;
}
//@{
/**
* Set the x, y, z and w components of a 3D vector in homogeneous coodinates.
*/
void Set(const T& x, const T& y, const T& z, const T& w)
{
this->Data[0] = x;
this->Data[1] = y;
this->Data[2] = z;
this->Data[3] = w;
}
//@}
/**
* Set the x component of the vector, i.e. element 0.
*/
void SetX(const T& x) { this->Data[0] = x; }
/**
* Get the x component of the vector, i.e. element 0.
*/
const T& GetX() const { return this->Data[0]; }
/**
* Set the y component of the vector, i.e. element 1.
*/
void SetY(const T& y) { this->Data[1] = y; }
/**
* Get the y component of the vector, i.e. element 1.
*/
const T& GetY() const { return this->Data[1]; }
/**
* Set the z component of the vector, i.e. element 2.
*/
void SetZ(const T& z) { this->Data[2] = z; }
/**
* Get the z component of the vector, i.e. element 2.
*/
const T& GetZ() const { return this->Data[2]; }
/**
* Set the w component of the vector, i.e. element 3.
*/
void SetW(const T& w) { this->Data[3] = w; }
/**
* Get the w component of the vector, i.e. element 3.
*/
const T& GetW() const { return this->Data[3]; }
};
//@}
/**
* Some inline functions for the derived types.
*/
......@@ -385,5 +466,15 @@ public:
vtkVector3Cross(vtkVector3d, double)
};
class vtkVector4d : public vtkVector4<double>
{
public:
using Superclass = vtkVector4<double>;
vtkVector4d() {}
vtkVector4d(double x, double y, double z, double w) :
vtkVector4<double>(x, y, z, w) {};
vtkVectorDerivedMacro(vtkVector4d, double, 4)
};
#endif // vtkVector_h
// VTK-HeaderTest-Exclude: vtkVector.h
......@@ -57,6 +57,9 @@ public:
double *GetBounds() VTK_SIZEHINT(6) override;
void GetBounds(double bounds[6]) override
{ this->vtkAbstractMapper3D::GetBounds(bounds); };
virtual double* GetBounds(const int /* port */)
{ return this->GetBounds(); };
//@}
//@{
......
......@@ -66,8 +66,8 @@ public:
/**
* Set/Get the volume property.
*/
void SetProperty(vtkVolumeProperty *property);
vtkVolumeProperty *GetProperty();
virtual void SetProperty(vtkVolumeProperty *property);
virtual vtkVolumeProperty *GetProperty();
//@}
/**
......@@ -202,10 +202,13 @@ public:
void UpdateScalarOpacityforSampleSize(vtkRenderer *ren,
float sample_distance);
/// Used by vtkHardwareSelector to determine if the prop supports hardware
/// selection.
/// @warning INTERNAL METHOD - NOT INTENDED FOR GENERAL USE
/// DO NOT USE THIS METHOD OUTSIDE OF THE RENDERING PROCESS
/**
* Used by vtkHardwareSelector to determine if the prop supports hardware
* selection.
*
* @warning INTERNAL METHOD - NOT INTENDED FOR GENERAL USE
* DO NOT USE THIS METHOD OUTSIDE OF THE RENDERING PROCESS
*/
bool GetSupportsSelection() override
{ return true; }
......
......@@ -395,32 +395,33 @@ bool vtkDualDepthPeelingPass::PreReplaceVolumetricShaderValues(
{
const std::string rayInit =
" // Transform zStart and zEnd to texture_coordinates\n"
" mat4 NDCToTextureCoords = ip_inverseTextureDataAdjusted * in_inverseVolumeMatrix *\n"
" mat4 NDCToTextureCoords = ip_inverseTextureDataAdjusted * in_inverseVolumeMatrix[0] *\n"
" in_inverseModelViewMatrix * in_inverseProjectionMatrix;\n"
" \n"
" // Start point\n"
" vec4 startPoint = WindowToNDC(gl_FragCoord.x, gl_FragCoord.y, zStart);\n"
" startPoint = NDCToTextureCoords * startPoint;\n"
" startPoint /= startPoint.w;\n"
" \n"
" // startPoint could be located outside of the bounding box (bbox), this\n"
" // is the case in:\n"
" // 1. PeelVolumesOutside: Areas external to any geometry.\n"
" // 2. PeelVolumetricGeometry: Areas where the volume is contained within\n"
" // translucent geometry but the containing geometry lies outside of the bbox\n"
" // (startPoint is either in-front or behind the bbox depending on the viewpoint).\n"
"\n"
" // Given that startPoint could be located either in-front, inside or behind the\n"
" // bbox (the ray exit is unknown hence it is not possible to use clamp() directly),\n"
" // the clamp is divided in these three zones:\n"
" // a. In-front: clamp to ip_textureCoords (bbox's texture coord).\n"
" // b. Inside: use startPoint directly as it is peeling within the bbox.\n"
" // c. Behind: discard by returning vec4(0.f).\n"
// startPoint could be located outside of the bounding box (bbox), this
// is the case in:
// 1. PeelVolumesOutside: Areas external to any geometry.
// 2. PeelVolumetricGeometry: Areas where the volume is contained within
// translucent geometry but the containing geometry lies outside of the bbox
// (startPoint is either in-front or behind the bbox depending on the viewpoint).
//
// Given that startPoint could be located either in-front, inside or behind the\n"
// bbox (the ray exit is unknown hence it is not possible to use clamp() directly),\n"
// the clamp is divided in these three zones:\n"
// a. In-front: clamp to ip_textureCoords (bbox's texture coord).\n"
// b. Inside: use startPoint directly as it is peeling within the bbox.\n"
// c. Behind: discard by returning vec4(0.f).\n"
"\n"
" // Initialize g_dataPos as if startPoint lies Inside (b.)\n"
" g_dataPos = startPoint.xyz;\n"
" bool isInsideBBox = !(any(greaterThan(startPoint.xyz, in_texMax)) ||\n"
" any(lessThan(startPoint.xyz, in_texMin)));\n"
" bool isInsideBBox = !(any(greaterThan(startPoint.xyz, in_texMax[0])) ||\n"
" any(lessThan(startPoint.xyz, in_texMin[0])));\n"
" if (!isInsideBBox)\n"
" {\n"
" vec3 distStartTexCoord = ip_textureCoords.xyz - startPoint.xyz;\n"
......@@ -602,23 +603,27 @@ bool vtkDualDepthPeelingPass::PreReplaceVolumetricShaderValues(
" bool onlyBack = frontStartDepth == backStartDepth &&\n"
" frontEndDepth == backEndDepth;\n"
"\n"
" // In the last peel, innerDepths may be (-1, -1) for most of the\n"
" // fragments. Casting a ray from [outerDepths.x, 1.0] would result\n"
" // in accumulating areas that have already been accounted for in\n"
" // former volume peels. In this case frontEndDepth should be the\n"
" // outer max instead. Because of this, the back castRay() is also\n"
" // skipped.\n"
// In the last peel, innerDepths may be (-1, -1) for most of the
// fragments. Casting a ray from [outerDepths.x, 1.0] would result
// in accumulating areas that have already been accounted for in
// former volume peels. In this case frontEndDepth should be the
// outer max instead. Because of this, the back castRay() is also
// skipped.
" bool noInnerDepths = innerDepths.x == -1.0;\n"
" if (noInnerDepths)\n"
" {\n"
" frontEndDepth = outerDepths.y;\n"
" }\n"
"\n"
" // Peel passes set -1 in pixels that contain only opaque geometry,\n"
" // so the opaque depth is fetched in order to z-composite volumes\n"
" // with opaque goemetry. To do this, the end point of front is clamped\n"
" // to opaque-depth and back ray-cast is skipped altogether since it\n"
" // would be covered by opaque geometry anyway.\n"
// Peel passes set -1 in pixels that contain only opaque geometry,
// so the opaque depth is fetched in order to z-composite volumes
// with opaque goemetry. To do this, the end point of front is clamped
// to opaque-depth and back ray-cast is skipped altogether since it
// would be covered by opaque geometry anyway.
" float oDepth = texture2D(opaqueDepthTex, pixelCoord * in_inverseWindowSize).x;\n"
" bool endBehindOpaque = frontEndDepth >= oDepth;\n"
" float clampedFrontEnd = frontEndDepth;\n"
......
......@@ -12,6 +12,7 @@ set(Module_SRCS
vtkFixedPointVolumeRayCastMIPHelper.cxx
vtkFixedPointVolumeRayCastMapper.cxx
vtkGPUVolumeRayCastMapper.cxx
vtkMultiVolume.cxx
vtkOSPRayVolumeInterface.cxx
vtkProjectedTetrahedraMapper.cxx
vtkRayCastImageDisplayHelper.cxx
......
......@@ -55,6 +55,10 @@ set (VolumeCxxTests
TestGPURayCastJitteringCustom.cxx
TestGPURayCastLargeColorTransferFunction.cxx
TestGPURayCastMapperSampleDistance.cxx
TestGPURayCastMultiVolumeAddRemove.cxx
TestGPURayCastMultiVolumeCellData.cxx
TestGPURayCastMultiVolumeOverlapping.cxx
TestGPURayCastMultiVolumeTransfer2D.cxx
TestGPURayCastPositionalLights.cxx
TestGPURayCastReleaseResources.cxx
TestGPURayCastRenderDepthToImage.cxx
......
/*=========================================================================
Program: Visualization Toolkit
Module: TestGPURayCastMultiVolumeAddRemove.cxx
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
All rights reserved.
See Copyright.txt or http://www.kitware.com/Copyright.htm 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.
=========================================================================*/
/**
* Tests adding and removing inputs to vtkMultiVolume and vtkGPUVolumeRCMapper.
*
*/
#include "vtkAxesActor.h"
#include "vtkCamera.h"
#include "vtkColorTransferFunction.h"
#include "vtkCommand.h"
#include "vtkConeSource.h"
#include "vtkGPUVolumeRayCastMapper.h"
#include "vtkOpenGLGPUVolumeRayCastMapper.h"
#include "vtkInteractorStyleTrackballCamera.h"
#include "vtkXMLImageDataReader.h"
#include "vtkVolume16Reader.h"
#include "vtkNew.h"
#include "vtkPiecewiseFunction.h"
#include "vtkRegressionTestImage.h"
#include "vtkRenderer.h"
#include "vtkRenderWindow.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkTestUtilities.h"
#include "vtkMultiVolume.h"
#include "vtkVolumeProperty.h"
#include "vtkImageResize.h"
#include "vtkImageResample.h"
#include "vtkImageData.h"
#include "vtkOutlineFilter.h"
#include "vtkPolyDataMapper.h"
#include "vtkAbstractMapper.h"
#include "vtkMath.h"
#include <chrono>
int TestGPURayCastMultiVolumeAddRemove(int argc, char* argv[])
{
// Load data
vtkNew<vtkVolume16Reader> reader;
reader->SetDataDimensions(64, 64);
reader->SetImageRange(1, 93);
reader->SetDataByteOrderToLittleEndian();
char* fname = vtkTestUtilities::ExpandDataFileName(argc, argv, "Data/headsq/quarter");
reader->SetFilePrefix(fname);
delete[] fname;
reader->SetDataSpacing(3.2, 3.2, 1.5);
vtkNew<vtkXMLImageDataReader> vaseSource;
const char* volumeFile = vtkTestUtilities::ExpandDataFileName(argc, argv, "Data/vase_1comp.vti");
vaseSource->SetFileName(volumeFile);
delete[] volumeFile;
vtkSmartPointer<vtkXMLImageDataReader> xmlReader =
vtkSmartPointer<vtkXMLImageDataReader>::New();
char* filename =
vtkTestUtilities::ExpandDataFileName(argc, argv, "Data/hncma-atlas.vti");
xmlReader->SetFileName(filename);
xmlReader->Update();
delete [] filename;
filename = nullptr;
// Volume 0 (upsampled headmr)
// ---------------------------
vtkNew<vtkImageResize> headmrSource;
headmrSource->SetInputConnection(reader->GetOutputPort());
headmrSource->SetResizeMethodToOutputDimensions();
headmrSource->SetOutputDimensions(128, 128, 128);
headmrSource->Update();
vtkNew<vtkColorTransferFunction> ctf;
ctf->AddRGBPoint(0, 0.0, 0.0, 0.0);
ctf->AddRGBPoint(500, 1.0, 0.5, 0.3);
ctf->AddRGBPoint(1000, 1.0, 0.5, 0.3);
ctf->AddRGBPoint(1150, 1.0, 1.0, 0.9);
vtkNew<vtkPiecewiseFunction> pf;
pf->AddPoint(0, 0.00);
pf->AddPoint(500, 0.15);
pf->AddPoint(1000, 0.15);
pf->AddPoint(1150, 0.85);
vtkNew<vtkPiecewiseFunction> gf;
gf->AddPoint(0, 0.0);
gf->AddPoint(90, 0.1);
gf->AddPoint(100, 0.7);
vtkNew<vtkVolume> vol;
vol->GetProperty()->SetScalarOpacity(pf);
vol->GetProperty()->SetColor(ctf);
vol->GetProperty()->SetGradientOpacity(gf);
vol->GetProperty()->SetInterpolationType(VTK_LINEAR_INTERPOLATION);
// Volume 1 (vase)
// -----------------------------
vtkNew<vtkColorTransferFunction> ctf1;
ctf1->AddRGBPoint(0, 0.0, 0.0, 0.0);
ctf1->AddRGBPoint(500, 0.1, 1.0, 0.3);
ctf1->AddRGBPoint(1000, 0.1, 1.0, 0.3);
ctf1->AddRGBPoint(1150, 1.0, 1.0, 0.9);
vtkNew<vtkPiecewiseFunction> pf1;
pf1->AddPoint(0, 0.0);
pf1->AddPoint(500, 1.0);
vtkNew<vtkPiecewiseFunction> gf1;
gf1->AddPoint(0, 0.0);
gf1->AddPoint(550, 1.0);
vtkNew<vtkVolume> vol1;
vol1->GetProperty()->SetScalarOpacity(pf1);
vol1->GetProperty()->SetColor(ctf1);
vol1->GetProperty()->SetGradientOpacity(gf1);
vol1->GetProperty()->SetInterpolationType(VTK_LINEAR_INTERPOLATION);
vol1->RotateX(-55.);
vol1->SetPosition(80., 50., 130.);
// Volume 2 (brain)
// -----------------------------
vtkNew<vtkPiecewiseFunction> pf2;
pf1->AddPoint(0, 0.0);
pf1->AddPoint(5022, 0.09);
vtkNew<vtkColorTransferFunction> ctf2;
ctf2->AddRGBPoint(0, 1.0, 0.3, 0.2);
ctf2->AddRGBPoint(2511, 0.3, 0.2, 0.9);
ctf2->AddRGBPoint(5022, 0.5, 0.6, 1.0);
vtkNew<vtkPiecewiseFunction> gf2;
gf2->AddPoint(0, 0.0);
gf2->AddPoint(550, 0.5);
vtkNew<vtkVolume> vol2;
vol2->GetProperty()->SetScalarOpacity(pf2);
vol2->GetProperty()->SetColor(ctf2);
vol2->GetProperty()->SetGradientOpacity(gf2);
vol2->GetProperty()->SetInterpolationType(VTK_LINEAR_INTERPOLATION);
vol2->SetScale(0.8, 0.8, 0.8);
vol2->SetPosition(210., 200., -90.);
vol2->RotateX(90.);
vol2->RotateY(-95.);
vol2->RotateZ(-5.);
// Rendering context
vtkNew<vtkRenderWindow> renWin;
renWin->SetSize(512, 512);
renWin->SetMultiSamples(0);
vtkNew<vtkRenderer> ren;
renWin->AddRenderer(ren);
ren->SetBackground(0.0, 0.0, 0.0);
vtkNew<vtkRenderWindowInteractor> iren;
iren->SetRenderWindow(renWin);
vtkNew<vtkInteractorStyleTrackballCamera> style;
iren->SetInteractorStyle(style);
auto cam = ren->GetActiveCamera();
cam->SetFocalPoint(41.9596, -17.9662, 78.5903);
cam->SetPosition(373.891, 619.954, -53.5932);
cam->SetViewUp(-0.0358384, -0.184856, -0.982112);
renWin->Render();
// Multi volume instance
// ---------------------
vtkNew<vtkMultiVolume> overlappingVol;
vtkNew<vtkGPUVolumeRayCastMapper> mapper;
mapper->SetUseJittering(0);
overlappingVol->SetMapper(mapper);
// Paramters that are global to all of the inputs are currently
// set through the vtkVolumeProperty corresponding to the required
// input port (port 0)
vol->GetProperty()->SetInterpolationType(VTK_LINEAR_INTERPOLATION);
mapper->SetInputConnection(4, xmlReader->GetOutputPort());
overlappingVol->SetVolume(vol2, 4);
mapper->SetInputConnection(0, headmrSource->GetOutputPort());
overlappingVol->SetVolume(vol, 0);
mapper->SetInputConnection(2, vaseSource->GetOutputPort());
overlappingVol->SetVolume(vol1, 2);
ren->AddVolume(overlappingVol);
renWin->Render();
// Remove / add
mapper->RemoveInputConnection(4, 0);
overlappingVol->RemoveVolume(4);
renWin->Render();
mapper->RemoveInputConnection(2, 0);
overlappingVol->RemoveVolume(2);
renWin->Render();
mapper->SetInputConnection(4, xmlReader->GetOutputPort());
overlappingVol->SetVolume(vol2, 4);
renWin->Render();
int retVal = vtkTesting::Test(argc, argv, renWin, 90);
if (retVal == vtkRegressionTester::DO_INTERACTOR)
{
iren->Start();
}
return !((retVal == vtkTesting::PASSED) ||
(retVal == vtkTesting::DO_INTERACTOR));
}
/*=========================================================================
Program: Visualization Toolkit
Module: TestGPURayCastMultiVolumeCellData.cxx
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
All rights reserved.
See Copyright.txt or http://www.kitware.com/Copyright.htm 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.
=========================================================================*/
/**
* Sets two inputs in vtkGPUVolumeRayCastMapper and uses a vtkMultiVolume
* instance to render the two inputs simultaneously (one point-data and one
* cell-data). Each vtkVolume contains independent transfer functions (one
* a set of 1D Tfs and the other a 2D Tf).
*/
#include "vtkAxesActor.h"
#include "vtkCamera.h"
#include "vtkPointDataToCellData.h"
#include "vtkColorTransferFunction.h"
#include "vtkCommand.h"
#include "vtkConeSource.h"
#include "vtkGPUVolumeRayCastMapper.h"
#include "vtkOpenGLGPUVolumeRayCastMapper.h"
#include "vtkInteractorStyleTrackballCamera.h"
#include "vtkXMLImageDataReader.h"
#include "vtkVolume16Reader.h"
#include "vtkNew.h"
#include "vtkNrrdReader.h"
#include "vtkPiecewiseFunction.h"
#include "vtkPNGReader.h"
#include "vtkPointData.h"
#include "vtkRegressionTestImage.h"
#include "vtkRenderer.h"
#include "vtkRenderWindow.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkTestUtilities.h"
#include "vtkMultiVolume.h"
#include "vtkVolumeProperty.h"
#include "vtkImageResize.h"
#include "vtkImageResample.h"
#include "vtkImageData.h"
#include "vtkOutlineFilter.h"
#include "vtkPolyDataMapper.h"
#include "vtkAbstractMapper.h"
#include "vtkMath.h"
#include <chrono>
namespace {
vtkSmartPointer<vtkImageData> ConvertImageToFloat(vtkDataObject* image)
{
auto imageIn = vtkImageData::SafeDownCast(image);
vtkSmartPointer<vtkImageData> imageOut = vtkSmartPointer<vtkImageData>::New();
imageOut->SetDimensions(imageIn->GetDimensions());
imageOut->AllocateScalars(VTK_FLOAT, 4);
auto arrayIn = imageIn->GetPointData()->GetScalars();
auto arrayOut = imageOut->GetPointData()->GetScalars();
const auto numTuples = arrayOut->GetNumberOfTuples();
for (vtkIdType i = 0; i < numTuples; i++)
{
double* value = arrayIn->GetTuple(i);
double valuef[4];
valuef[0] = value[0] / 255.;
valuef[1] = value[1] / 255.;
valuef[2] = value[2] / 255.;
valuef[3] = value[3] / 255.;
arrayOut->SetTuple(i, valuef);
}
//return std::move(imageOut);
return imageOut;
};
}
////////////////////////////////////////////////////////////////////////////////
int TestGPURayCastMultiVolumeCellData(int argc, char* argv[])
{
// Load data
vtkNew<vtkVolume16Reader> headReader;
headReader->SetDataDimensions(64, 64);
headReader->SetImageRange(1, 93);
headReader->SetDataByteOrderToLittleEndian();
char* fname = vtkTestUtilities::ExpandDataFileName(argc, argv,
"Data/headsq/quarter");
headReader->SetFilePrefix(fname);
headReader->SetDataSpacing(3.2, 3.2, 1.5);
delete[] fname;
fname = vtkTestUtilities::ExpandDataFileName(argc, argv,
"Data/tooth.nhdr");
vtkNew<vtkNrrdReader> toothReader;
toothReader->SetFileName(fname);
delete[] fname;
vtkNew<vtkPNGReader> reader2dtf;
fname = vtkTestUtilities::ExpandDataFileName(argc,
argv, "Data/tooth_2dtransf.png");
reader2dtf->SetFileName(fname);
reader2dtf->Update();
delete[] fname;
vtkNew<vtkAxesActor> axis;
axis->SetTotalLength(100., 100., 100.);
axis->SetNormalizedTipLength(0.1, 0.1, 0.1);
axis->SetNormalizedShaftLength(1., 1., 1.);
axis->AxisLabelsOff();
axis->SetConeRadius(0.5);
// Volume 0 (upsampled headmr)
// ---------------------------
// Transform the head dataset to cells
vtkNew<vtkImageResize> headmrSource;
headmrSource->SetInputConnection(headReader->GetOutputPort());
headmrSource->SetResizeMethodToOutputDimensions();
headmrSource->SetOutputDimensions(128, 128, 128);
vtkNew<vtkPointDataToCellData> pointsToCells;
pointsToCells->SetInputConnection(headmrSource->GetOutputPort());
pointsToCells->Update();
vtkNew<vtkColorTransferFunction> ctf;
ctf->AddRGBPoint(0, 0.0, 0.0, 0.0);
ctf->AddRGBPoint(500, 0.1, 0.6, 0.3);
ctf->AddRGBPoint(1000, 0.1, 0.6, 0.3);
ctf->AddRGBPoint(1150, 1.0, 1.0, 0.9);
vtkNew<vtkPiecewiseFunction> pf;
pf->AddPoint(0, 0.00);
pf->AddPoint(500, 0.15);
pf->AddPoint(1000, 0.15);
pf->AddPoint(1150, 0.85);
vtkNew<vtkPiecewiseFunction> gf;
gf->AddPoint(0, 0.0);
gf->AddPoint(90, 0.07);
gf->AddPoint(100, 0.7);
vtkNew<vtkVolume> vol;
vol->GetProperty()->SetScalarOpacity(pf);
vol->GetProperty()->SetColor(ctf);
vol->GetProperty()->SetGradientOpacity(gf);
vol->GetProperty()->SetInterpolationType(VTK_LINEAR_INTERPOLATION);
//vol->GetProperty()->ShadeOn();
// Volume 1 (tooth)
// -----------------------------
vtkNew<vtkVolume> vol1;
auto tf2d = ConvertImageToFloat(reader2dtf->GetOutputDataObject(0));
vol1->GetProperty()->SetTransferFunction2D(tf2d);
vol1->GetProperty()->SetInterpolationType(VTK_LINEAR_INTERPOLATION);
vol1->RotateX(180.);
vol1->RotateZ(90.);
vol1->SetScale(1.8, 1.8, 1.8);
vol1->SetPosition(175., 190., 210.);
// Multi volume instance
// ---------------------
// Create an overlapping volume prop (add specific properties to each
// entity).
vtkNew<vtkMultiVolume> overlappingVol;
vtkNew<vtkGPUVolumeRayCastMapper> mapper;
overlappingVol->SetMapper(mapper);
mapper->SetInputConnection(0, pointsToCells->GetOutputPort());
overlappingVol->SetVolume(vol, 0);
mapper->SetInputConnection(3, toothReader->GetOutputPort());
overlappingVol->SetVolume(vol1, 3);
mapper->SetUseJittering(1);
// Rendering context
vtkNew<vtkRenderWindow> renWin;
renWin->SetSize(800, 400);
renWin->SetMultiSamples(0);
// Outside renderer (left)
vtkNew<vtkRenderer> ren;
renWin->AddRenderer(ren);
ren->SetBackground(1.0, 1.0, 1.0);
ren->SetViewport(0.0, 0.0, 0.5, 1.0);
ren->AddActor(axis);
ren->AddVolume(overlappingVol);
vtkNew<vtkRenderWindowInteractor> iren;
iren->SetRenderWindow(renWin);
vtkNew<vtkInteractorStyleTrackballCamera> style;
iren->SetInteractorStyle(style);