学习整理黑山老妖博客《VTK笔记-切割分割-三维曲线或几何切割数据，通过阈值填充获得切割数据的结果。URL:https://blog.csdn.net/liushao1031177/article/details/118419221
我测试的图形如下：
整理后的源代码如下：
#pragma once #include <vtkNIFTIImageReader.h> #include <vtkInformation.h> #include <vtkVertexGlyphFilter.h> #include <vtkSelectEnclosedPoints.h> #include <vtkFloatArray.h> #include <vtkShortArray.h> #include <vtkActor.h> #include <vtkDICOMImageReader.h> #include <vtkPlane.h> #include <vtkPolygon.h> #include <vtkPolyDataMapper.h> #include <vtkGPUVolumeRayCastMapper.h> #include <vtkPolyDataToImageStencil.h> #include <vtkImageData.h> #include <vtkImageStencilToImage.h> #include <vtkLineSource.h> #include <vtkAutoInit.h> #include <vtkDataSet.h> #include <vtkPointData.h> #include <vtkVolumeProperty.h> #include <vtkPiecewiseFunction.h> #include <vtkVolume.h> #include <vtkColorTransferFunction.h> #include <vtkRenderer.h> #include <vtkRenderWindow.h> #include <vtkRenderWindowInteractor.h> #include <vtkPlaneSource.h> #include <vtkCamera.h> #include <vtkProperty.h> VTK_MODULE_INIT(vtkRenderingOpenGL2); VTK_MODULE_INIT(vtkRenderingVolumeOpenGL2); VTK_MODULE_INIT(vtkRenderingFreeType); VTK_MODULE_INIT(vtkInteractionStyle); double* spacing; double* origin; int* dm; class ClipVolumeByMask { public: static void Test() { vtkNew<vtkRenderer> ren; vtkNew<vtkDICOMImageReader> reader; reader->SetDirectoryName("D:\\VTK-DATA\\Dicom-Samples\\xuebaolin"); reader->Update(); vtkImageData* pImageData = reader->GetOutput(); dm = pImageData->GetDimensions(); spacing = reader->GetDataSpacing(); origin = reader->GetDataOrigin(); int extent[6] = { 0, dm[0], 0, dm[1], 0, dm[2] - 1 }; double camera_position[3] = { origin[0] spacing[0] * dm[0] / 2, origin[1] - spacing[1] * dm[1] * 2 , origin[2] spacing[2] * dm[2] / 2 }; double top_left[5][3] = { { camera_position[0] - 30,camera_position[1],camera_position[2] - 20 }, { camera_position[0] 30,camera_position[1],camera_position[2] - 20 }, { camera_position[0] 30,camera_position[1],camera_position[2] 20 }, { camera_position[0] - 30,camera_position[1],camera_position[2] 20 }, { camera_position[0], camera_position[1],camera_position[2]} }; double normal_camera[3] = { 0.0000000000000000, -1.0000000000000000, 0.0000000000000000 }; double another_polygon[5][3] = { 0 }; Get_Another_Plane_Polygon(normal_camera, top_left, another_polygon); #pragma region 产生五棱柱体 vtkSmartPointer<vtkPolyData> cylinderCutterBody = GenerateFivePrismatic(top_left, another_polygon); vtkSmartPointer<vtkPolyDataMapper> cutterMapper = vtkSmartPointer<vtkPolyDataMapper>::New(); cutterMapper->SetInputData(cylinderCutterBody); vtkSmartPointer<vtkActor> cutterActor = vtkSmartPointer<vtkActor>::New(); cutterActor->SetMapper(cutterMapper); #pragma endregion #pragma region 裁剪 vtkNew<vtkImageData> whiteImage; whiteImage->SetSpacing(spacing); whiteImage->SetExtent(extent); whiteImage->SetOrigin(origin); whiteImage->AllocateScalars(VTK_UNSIGNED_CHAR, 1); vtkNew<vtkPolyDataToImageStencil> polyToStencil; // <-- polyToStencil->SetInputData(cylinderCutterBody); polyToStencil->SetOutputOrigin(origin); polyToStencil->SetOutputSpacing(spacing); polyToStencil->SetOutputWholeExtent(whiteImage->GetExtent()); polyToStencil->Update(); const unsigned char inval = 0; //纯黑 const unsigned char outval = 255; //纯白 vtkNew<vtkImageStencilToImage> stencilToImage; // <-- stencilToImage->SetInputConnection(polyToStencil->GetOutputPort()); stencilToImage->SetInsideValue(inval); // <-- stencilToImage->SetOutsideValue(outval); // <-- stencilToImage->SetOutputScalarType(VTK_UNSIGNED_CHAR); stencilToImage->Update(); #pragma endregion // 在两个平面上画五个点，前五边形 vtkSmartPointer<vtkActor> pointsActor = DrawPolyFivePoint(camera_position, top_left, another_polygon); vtkSmartPointer<vtkActor> linesActor = DrawPolyFiveLine(camera_position, top_left, another_polygon); // 绘制视平面 vtkSmartPointer<vtkActor> planeActor = DrawViewPlane(normal_camera, camera_position); vtkSmartPointer<vtkActor> volumeActor = DrawOutline(); #pragma endregion #pragma region 使用vtkGPUVolumeRayCastMapper体渲染方式 vtkNew<vtkGPUVolumeRayCastMapper> gpuMapper; gpuMapper->SetInputData(pImageData); gpuMapper->SetMaskTypeToBinary(); gpuMapper->SetMaskInput(stenclToImage->GetOutput());
vtkNew<vtkVolumeProperty> volumeProperty;
volumeProperty->SetInterpolationTypeToLinear();
volumeProperty->SetAmbient(0.4);
volumeProperty->SetDiffuse(1.6); //漫反射
volumeProperty->SetSpecular(0.2); //镜面反射
volumeProperty->ShadeOn(); //打开或者关闭阴影测试
//设置不透明度
vtkNew<vtkPiecewiseFunction> compositeOpacity;
compositeOpacity->AddPoint(220, 0.00);
compositeOpacity->AddPoint(255, 1.0);
volumeProperty->SetScalarOpacity(compositeOpacity); //设置不透明度传输函数
vtkNew<vtkColorTransferFunction> color;
color->AddRGBPoint(0.000, 0.00, 0.00, 0.00);
color->AddRGBPoint(60, 1.00, 0.52, 0.30);
color->AddRGBPoint(150, 1.00, 1.00, 1.00);
color->AddRGBPoint(200, .5, .3, .2);
volumeProperty->SetColor(color);
vtkNew<vtkVolume> volume;
volume->SetMapper(gpuMapper);
volume->SetProperty(volumeProperty);
#pragma endregion
#pragma region 渲染管线
ren->SetBackground(0.3, 0.4, 0.3);
ren->AddVolume(volume);
ren->GetActiveCamera()->SetViewUp(0, 0, 1);
ren->GetActiveCamera()->SetFocalPoint(0, 0, 0);
ren->GetActiveCamera()->SetPosition(0, -1, 0);
ren->ResetCamera();
ren->AddActor(planeActor);
ren->AddActor(linesActor);
ren->AddActor(pointsActor);
ren->AddActor(cutterActor);
ren->AddActor(volumeActor);
planeActor->SetVisibility(false);
cutterActor->SetVisibility(false);
vtkNew<vtkRenderWindow> rw;
rw->AddRenderer(ren);
rw->SetSize(640, 480);
rw->Render();
rw->SetWindowName("VolumeRendering PipeLine");
vtkNew<vtkRenderWindowInteractor> rwi;
rwi->SetRenderWindow(rw);
rw->Render();
rwi->Start();
#pragma endregion
}
private:
static void Get_Another_Plane_Polygon(double normal[3], double polygon[5][3], double another_polygon[5][3])
{
double max_distance = Cal_Max_Distance(normal, polygon[0]);
for (int i = 0; i < 5; i++) {
another_polygon[i][0] = polygon[i][0] - max_distance * normal[0];
another_polygon[i][1] = polygon[i][1] - max_distance * normal[1];
another_polygon[i][2] = polygon[i][2] - max_distance * normal[2];
}
}
static double Cal_Max_Distance(double normal[3], double origin[3])
{
vtkNew<vtkPlane> plane;
plane->SetNormal(normal);
plane->SetOrigin(origin);
double point_000[3] = { 0,0,0 };
double dis1 = plane->DistanceToPlane(point_000);
double point_100[3] = { spacing[0] * dm[0],0,0 };
double dis2 = plane->DistanceToPlane(point_100);
double point_010[3] = { 0,spacing[0] * dm[0],0 };
double dis3 = plane->DistanceToPlane(point_010);
double point_110[3] = { spacing[0] * dm[0],spacing[0] * dm[0],0 };
double dis4 = plane->DistanceToPlane(point_110);
double point_001[3] = { 0,0,(dm[2] - 1) };
double dis5 = plane->DistanceToPlane(point_001);
double point_101[3] = { spacing[0] * dm[0],0,(dm[2] - 1) };
double dis6 = plane->DistanceToPlane(point_101);
double point_011[3] = { 0,spacing[0] * dm[0],(dm[2] - 1) };
double dis7 = plane->DistanceToPlane(point_011);
double point_111[3] = { spacing[0] * dm[0],spacing[0] * dm[0],(dm[2] - 1) };
double dis8 = plane->DistanceToPlane(point_111);
return dis8;
}
static vtkSmartPointer<vtkActor> DrawViewPlane(double normal[3], double center[3])
{
vtkSmartPointer<vtkPlaneSource> plane = vtkSmartPointer<vtkPlaneSource>::New();
plane->SetNormal(normal);
plane->SetCenter(center);
plane->SetOrigin(center[0] - 100, center[1], center[2] - 100);
plane->SetPoint1(center[0] + 100, center[1], center[2] - 100);
plane->SetPoint2(center[0] - 100, center[1], center[2] + 100);
vtkSmartPointer<vtkPolyDataMapper> planeMapper = vtkSmartPointer<vtkPolyDataMapper>::New();
planeMapper->SetInputConnection(plane->GetOutputPort());
vtkSmartPointer<vtkActor> planeActor = vtkSmartPointer<vtkActor>::New();
planeActor->SetMapper(planeMapper);
return planeActor;
}
static vtkSmartPointer<vtkActor> DrawPolyFivePoint(double center[3], double top_left[4][3], double another_polygon[4][3])
{
vtkSmartPointer<vtkPoints> points1 = vtkSmartPointer<vtkPoints>::New();
points1->InsertNextPoint(center);
points1->InsertNextPoint(top_left[0]);
points1->InsertNextPoint(top_left[1]);
points1->InsertNextPoint(top_left[2]);
points1->InsertNextPoint(top_left[3]);
points1->InsertNextPoint(another_polygon[0]);
points1->InsertNextPoint(another_polygon[1]);
points1->InsertNextPoint(another_polygon[2]);
points1->InsertNextPoint(another_polygon[3]);
vtkSmartPointer<vtkPolyData> pointsPolydata = vtkSmartPointer<vtkPolyData>::New();
pointsPolydata->SetPoints(points1);
vtkSmartPointer<vtkVertexGlyphFilter> vertexGlyphFilter = vtkSmartPointer<vtkVertexGlyphFilter>::New();
vertexGlyphFilter->AddInputData(pointsPolydata);
vertexGlyphFilter->Update();
vtkSmartPointer<vtkPolyDataMapper> pointsMapper = vtkSmartPointer<vtkPolyDataMapper>::New();
pointsMapper->SetInputConnection(vertexGlyphFilter->GetOutputPort());
vtkSmartPointer<vtkActor> pointsActor = vtkSmartPointer<vtkActor>::New();
pointsActor->SetMapper(pointsMapper);
pointsActor->GetProperty()->SetPointSize(10);//定义点的尺寸大小，这样点才能在画布上显示出来
pointsActor->GetProperty()->SetColor(0.0, 0.0, 1);
return pointsActor;
}
static vtkSmartPointer<vtkActor> DrawPolyFiveLine(double center[3], double top_left[4][3], double another_polygon[4][3])
{
vtkSmartPointer<vtkPoints> points1 = vtkSmartPointer<vtkPoints>::New();
points1->InsertNextPoint(center);
points1->InsertNextPoint(top_left[0]);
points1->InsertNextPoint(top_left[1]);
points1->InsertNextPoint(top_left[2]);
points1->InsertNextPoint(top_left[3]);
vtkSmartPointer<vtkCellArray> line_polys = vtkSmartPointer<vtkCellArray>::New();
line_polys->InsertNextCell(6);
for (size_t i = 0; i < 5; i++)
{
line_polys->InsertCellPoint(i);
}
line_polys->InsertCellPoint(0);
vtkSmartPointer<vtkPolyData> geometry_line = vtkSmartPointer<vtkPolyData>::New();
geometry_line->SetPoints(points1);
geometry_line->SetLines(line_polys);
vtkSmartPointer<vtkPolyDataMapper> linesMapper = vtkSmartPointer<vtkPolyDataMapper>::New();
linesMapper->SetInputData(geometry_line);
vtkSmartPointer<vtkActor> linesActor = vtkSmartPointer<vtkActor>::New();
linesActor->SetMapper(linesMapper);
linesActor->GetProperty()->SetColor(255, 0, 0);
return linesActor;
}
static vtkSmartPointer<vtkPolyData> GenerateFivePrismatic(double top_left[4][3], double another_polygon[4][3])
{
vtkSmartPointer<vtkPolyData> geometry = vtkSmartPointer<vtkPolyData>::New();
vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New();
for (size_t i = 0; i < 5; i++)
{
points->InsertPoint(2 * i, top_left[i]);
points->InsertPoint(2 * i + 1, another_polygon[i]);
}
vtkSmartPointer<vtkCellArray> polys = vtkSmartPointer<vtkCellArray>::New();
vtkSmartPointer<vtkCellArray> strips = vtkSmartPointer<vtkCellArray>::New();
// 前平面
polys->InsertNextCell(5);
for (int i = 0; i < 5; i++) {
polys->InsertCellPoint(2 * i);
}
// 后平面
polys->InsertNextCell(5);
for (int i = 0; i < 5; i++) {
polys->InsertCellPoint(2 * i + 1);
}
// 中间柱面
strips->InsertNextCell(10 + 2);
for (int i = 0; i < 10; i++)
strips->InsertCellPoint(i);
strips->InsertCellPoint(0);
strips->InsertCellPoint(1);
geometry->SetPoints(points);
geometry->SetPolys(polys);
geometry->SetStrips(strips);
return geometry;
}
static vtkSmartPointer<vtkActor> DrawOutline()
{
vtkSmartPointer<vtkPoints> points_volume = vtkSmartPointer<vtkPoints>::New();
points_volume->InsertNextPoint(0, 0, 0);
points_volume->InsertNextPoint(spacing[0] * dm[0], 0, 0);
points_volume->InsertNextPoint(spacing[0] * dm[0], spacing[0] * dm[0], 0);
points_volume->InsertNextPoint(0, spacing[0] * dm[0], 0);
points_volume->InsertNextPoint(0, 0, (dm[2] - 1));
points_volume->InsertNextPoint(spacing[0] * dm[0], 0, (dm[2] - 1));
points_volume->InsertNextPoint(spacing[0] * dm[0], spacing[0] * dm[0], (dm[2] - 1));
points_volume->InsertNextPoint(0, spacing[0] * dm[0], (dm[2] - 1));
vtkSmartPointer<vtkCellArray> line_volume = vtkSmartPointer<vtkCellArray>::New();
line_volume->InsertNextCell(5);
for (int i = 0; i < 4; i++)
line_volume->InsertCellPoint(i);
line_volume->InsertCellPoint(0);
line_volume->InsertNextCell(5);
for (int i = 4; i < 8; i++)
line_volume->InsertCellPoint(i);
line_volume->InsertCellPoint(4);
for (int i = 0; i < 4; i++)
{
line_volume->InsertNextCell(2);
line_volume->InsertCellPoint(0 + i);
line_volume->InsertCellPoint(4 + i);
}
vtkSmartPointer<vtkPolyData> geometry_volume = vtkSmartPointer<vtkPolyData>::New();
geometry_volume->SetPoints(points_volume);
geometry_volume->SetLines(line_volume);
vtkSmartPointer<vtkPolyDataMapper> volumeMapper = vtkSmartPointer<vtkPolyDataMapper>::New();
volumeMapper->SetInputData(geometry_volume);
vtkSmartPointer<vtkActor> volumeActor = vtkSmartPointer<vtkActor>::New();
volumeActor->SetMapper(volumeMapper);
return volumeActor;
}
};
int main(int, char* [])
{
ClipVolumeByMask::Test();
return 0;
}