visualize three-dimensional medical images

Development
1

Hello everyone,

I am using VTK to visualize three-dimensional medical images reconstructed from a series of CT slice images. I have attempted both volume rendering and surface rendering, but the results are not satisfactory. My ultimate goal is to visualize the entire brain in 3D. I suspect that the issue may be due to the limited number of slices. Could you please take a look and offer some advice?

Below is my code and some sample two-dimensional slice images.

Thank you!

image
image1902×1148 218 KB

#include <vtkSmartPointer.h>
#include <vtkImageData.h>
#include <vtkStructuredPoints.h>
#include <vtkStructuredPointsReader.h>
#include <vtkColorTransferFunction.h>
#include <vtkPiecewiseFunction.h>
#include <vtkRenderer.h>
#include <vtkRenderWindow.h>
#include <vtkRenderWindowInteractor.h>
#include <vtkVolumeProperty.h>
#include <vtkFixedPointVolumeRayCastMapper.h>
#include <initializer_list>
#include <vtkPolyData.h>
#include <vtkPolyDataReader.h>
#include <vtkAutoInit.h>
#include “vtkNIFTIImageReader.h”
#include “vtkStringArray.h”
#include “vtkTIFFReader.h”
#include “vtkImageGaussianSmooth.h”
#include “vtkImageLuminance.h”
VTK_MODULE_INIT(vtkRenderingOpenGL2) // Initialize OpenGL2 rendering module
VTK_MODULE_INIT(vtkRenderingVolumeOpenGL2) // Initialize OpenGL2 volume rendering module
VTK_MODULE_INIT(vtkRenderingFreeType) // Initialize FreeType rendering module
VTK_MODULE_INIT(vtkInteractionStyle) // Initialize interaction style module

int main(int argc, char* argv)
{
vtkSmartPointer fileArray = vtkSmartPointer ::New();
char fileName[128];
for (int i = 1; i < 21; i++) {
sprintf_s(fileName, “C:\Users\Lenovo\Desktop\TCGA_CS_4941_19960909\TCGA_CS_4942_19970222\TCGA_CS_4942_19970222_%d.tif”, i);

	// "C:\\Users\\Lenovo\\Desktop\\TCGA_CS_4941_19960909\\TCGA_CS_4941_19960909\\TCGA_CS_4941_19960909_%d_mask.tif
	//"C:\\Users\\Lenovo\\Desktop\\TCGA_CS_4941_19960909\\TCGA_CS_4942_19970222\\TCGA_CS_4942_19970222_%d_mask.tif"
	fileArray->InsertNextValue(fileName);
}

vtkSmartPointer<vtkTIFFReader>reader = vtkSmartPointer<vtkTIFFReader>::New();
reader->SetFileNames(fileArray);
reader->Update();

vtkSmartPointer<vtkImageLuminance> luminanceFilter = vtkSmartPointer<vtkImageLuminance>::New();
luminanceFilter->SetInputData(reader->GetOutput());
luminanceFilter->Update();


// Gaussian smoothing
vtkSmartPointer<vtkImageGaussianSmooth> gaussianSmooth = vtkSmartPointer<vtkImageGaussianSmooth>::New();
gaussianSmooth->SetInputData(luminanceFilter->GetOutput());
gaussianSmooth->SetStandardDeviation(1.0); // Set the standard deviation to control the degree of smoothing
gaussianSmooth->Update();

// Create volume rendering mapper
vtkNew<vtkFixedPointVolumeRayCastMapper> volumeMapper;
volumeMapper->SetInputData(gaussianSmooth->GetOutput()); // Set the smoothed data

vtkNew<vtkVolumeProperty> volumeProperty;
volumeProperty->SetInterpolationTypeToLinear(); // Linear interpolation
volumeProperty->SetAmbient(5.0);  // Ambient light coefficient
volumeProperty->SetDiffuse(0.5);  // Diffuse reflection light coefficient
volumeProperty->SetSpecular(0.5); // Specular reflection light coefficient

vtkNew<vtkPiecewiseFunction> compositeOpacity;
compositeOpacity->AddPoint(70, 0.00);
compositeOpacity->AddPoint(90, 0.20);
compositeOpacity->AddPoint(110, 0.40);
compositeOpacity->AddPoint(130, 0.60);
compositeOpacity->AddPoint(150, 0.80);
compositeOpacity->AddPoint(170, 1.00);
volumeProperty->SetScalarOpacity(compositeOpacity); // Set the opacity transfer function

vtkNew<vtkPiecewiseFunction> volumeGradientOpacity;
volumeGradientOpacity->AddPoint(10, 0.5);
volumeGradientOpacity->AddPoint(100, 1.0);
//volumeProperty->SetGradientOpacity(volumeGradientOpacity);// Set gradient opacity for comparison

vtkNew<vtkColorTransferFunction> color;
color->AddRGBPoint(0.000, 0.00, 0.00, 0.00);
color->AddRGBPoint(64.00, 1.80, 1.80, 0.30);
color->AddRGBPoint(190.0, 1.00, 1.00, 1.00);
color->AddRGBPoint(220.0, 0.20, 0.20, 0.20);
volumeProperty->SetColor(color);

vtkNew<vtkVolume> volume;
volume->SetMapper(volumeMapper);
volume->SetProperty(volumeProperty);


vtkNew<vtkRenderer> renderer;
renderer->SetBackground(0.0, 0.0, 0.0);
renderer->AddVolume(volume);



vtkNew<vtkRenderWindow> renWin;
renWin->AddRenderer(renderer);
renWin->SetSize(640, 480);
renWin->SetWindowName("raycasting_algorithm");

vtkNew<vtkRenderWindowInteractor> iren;
iren->SetRenderWindow(renWin);
renWin->Render();
iren->Start();

}