如何在 vtk 中的 3D 表面渲染输出中获取 2D dicom 图像切片的位置
How to get the position of 2D dicom image slice in 3D surface rendered output in vtk
我正在做一个 VTK 程序,当我输入 2D DICOM 患者图像位置 (请参考给定的图像以便更好地理解),我需要在 3D 表面渲染输出中获取特定切片。
对于体积渲染的 3D 图像,这可以通过使用这些函数来实现,即 vtkImageData、vtkImageMapToColors、vtkImage Actor。
我的问题是如何在渲染的表面上做到这一点 output.do 任何人都知道这个概念。如果有人知道请回答。如果我的问题不正确或难以理解,请分享您的意见。
为了便于理解,我展示了一张示例图片
考虑一下我的 3d 输出将如下图所示
和当我在文本框中输入图像位置(患者)并单击确定按钮时,对应的切片应显示在 3d 图像中,如下图所示
如果我的问题无法理解,请告知
我被困在这里我什至不知道我的代码是 right.Here 是不是我的代码
ExtractVOI->SetInput(reader1->GetOutput());//VOI extractor
ExtractVOI->SetVOI(1,0,0,0,1,0);//i have given the Image Orientation(patient) as the SetVOI value
////====CREATE LookUpTable
tableax1->SetTableRange(0.0, 4096.0);
tableax1->SetValueRange(0.0, 1.0);
tableax1->SetSaturationRange(0.0, 0.0);
tableax1->SetRampToSCurve();
tableax1->SetAlphaRange(0.0, 0.08);
tableax1->Build();
//====CREATE TEXTURE
planesourceax1->SetXResolution(1);
planesourceax1->SetYResolution(1);
planesourceax1->SetOrigin(0,0,0);
planesourceax1->SetPoint1(xg , yg,zg);//i have given the value of Image Position(patient) that is taken from a textbox ,as the points
planesourceax1->Update();
vtkSmartPointer<vtkPolyDataMapper> mapax1 = vtkSmartPointer<vtkPolyDataMapper>::New();
mapax1->SetInputConnection(planesourceax1->GetOutputPort());
mapax1->UpdateWholeExtent();
textureax1->SetInputConnection(ExtractVOI->GetOutputPort());
textureax1->InterpolateOn();
textureax1->SetLookupTable(tableax1);
textureax1->UpdateWholeExtent();
//===PASS TO ACTOR
actorax1->SetMapper(mapax1);
actorax1->GetMapper()->SetResolveCoincidentTopologyToPolygonOffset();
actorax1->GetMapper()->SetResolveCoincidentTopologyPolygonOffsetParameters(0.1, -1.0);
actorax1->SetTexture(textureax1);
renderer->AddActor(actorax1);
renderWindow->Render();
但我没有得到输出
我也试过了:
static double axialElements[16] = {
1, 0, 0, 0,
0, 1, 0, 0,
0, 0, 1, 0,
0, 0, 0, 1 };
resliceax1->SetInputConnection(reader->GetOutputPort());
resliceax1->SetOutputDimensionality(2);
vtkSmartPointer<vtkMatrix4x4> reslicematrixax1 = vtkSmartPointer<vtkMatrix4x4>::New();
reslicematrixax1->DeepCopy(axialElements);
resliceax1->SetResliceAxes(reslicematrixax1);
resliceax1->SetResliceAxesOrigin(0.0, 0.0, 0.0);
resliceax1->Update();
extractaxpos1->RemoveAllInputs();
extractaxpos1->SetInputConnection(resliceax1->GetOutputPort());
////====CREATE LUT
tableax1->SetTableRange(0.0, 4096.0);
tableax1->SetValueRange(0.0, 1.0);
tableax1->SetSaturationRange(0.0, 0.0);
tableax1->SetRampToSCurve();
tableax1->SetAlphaRange(0.0, 0.08);
tableax1->Build();
//====CREATE TEXTURE
planesourceax1->SetXResolution(1);
planesourceax1->SetYResolution(1);
planesourceax1->SetOrigin(0,0,0);
planesourceax1->SetPoint1((xval/20 + xval/32),(yval/20 + yval/32),(zval/20 + zval/32));//this is where i put the values ad divided by its tag id(0020,0032)
//planesourceax1->SetPoint2(fBounds[0] , fBounds[3], fBounds[4]);
planesourceax1->Update();
vtkSmartPointer<vtkPolyDataMapper> mapax1 = vtkSmartPointer<vtkPolyDataMapper>::New();
mapax1->SetInputConnection(planesourceax1->GetOutputPort());
mapax1->UpdateWholeExtent();
textureax1->SetInputConnection(extractaxpos1->GetOutputPort());
textureax1->InterpolateOn();
textureax1->SetLookupTable(tableax1);
textureax1->UpdateWholeExtent();
//===PASS TO ACTOR
actorax1->SetMapper(mapax1);
actorax1->GetMapper()->SetResolveCoincidentTopologyToPolygonOffset();
actorax1->GetMapper()->SetResolveCoincidentTopologyPolygonOffsetParameters(0.1, -1.0);
actorax1->SetTexture(textureax1);
resliceax1->SetResliceAxesOrigin(0.0, 0.0, 0.0);
actorax1->SetPosition((xval/20 + xval/32),(yval/20 + yval/32),(zval/20 + zval/32));//I made the same changes here also
planesourceax1->SetOrigin(fBoundsUpdated[0], fBoundsUpdated[2], pDoc->fBounds[4]);
planesourceax1->SetPoint1(fBoundsUpdated[1] , fBoundsUpdated[2], pDoc->fBounds[4]);
planesourceax1->SetPoint2(fBoundsUpdated[0] , fBoundsUpdated[3], pDoc->fBounds[4]);
planesourceax1->Update();
但不是切割切片的位置is.It是切割不同的位置
如果你想知道什么元素决定切割方向
https://github.com/Kitware/VTK/blob/master/Examples/ImageProcessing/Cxx/ImageSlicing.cxx
请看这个link
static double axialElements[16] = {
1, 0, 0, 0,
0, 1, 0, 0,
0, 0, 1, 0,
0, 0, 0, 1 };
...
vtkSmartPointer<vtkMatrix4x4> resliceAxes = vtkSmartPointer<vtkMatrix4x4>::New();
resliceAxes->DeepCopy(axialElements);
// Set the point through which to slice
resliceAxes->SetElement(0, 3, 0);
resliceAxes->SetElement(1, 3, 0);
resliceAxes->SetElement(2, 3, 70);
矩阵 4x4 正在设置 重新切片方向
和 SetElement 决定 Reslice origin
此代码设置选项为Reslice XYimage start at Z = 70
vtkSmartPointer<vtkImageReslice> reslice = vtkSmartPointer<vtkImageReslice>::New();
reslice->SetInputData(m_vtkVolumeImageData);
reslice->SetOutputDimensionality(2);
reslice->SetResliceAxes(resliceAxes);
reslice->SetInterpolationModeToLinear();
reslice->Update();
此代码生成 Reslice 数据。
结果维度为 2D(平面)
如果你设置这个顺序,
将 Reslice 数据连接到 vtkImageMapToColors 以绘制 Reslice 图像。
最后,连接 Mapper 和 Actor 以使用 Volume 进行展示。
您需要设置 Actor 位置,因为 Reslice Data 可能没有位置信息
如果平面位置没有改变,使用 vtkDataSetMapper 和 vtkActor,而不是 vtkImageActor
或者你只是使用Widget,
我推荐 vtkImagePlaneWidget。
操作简单,功能强大。
希望对您有所帮助。
如果需要完整代码,请告诉我
void CMFCApplication1Dlg::setSliceImage()
{
int extent[6];
double spacing[3];
double origin[3];
//m_vtkVolumeImageData is vtkImageData(Save Slice Images)
//Same as dicomReader->GetOutput();
m_vtkVolumeImageData->GetExtent(extent);
m_vtkVolumeImageData->GetSpacing(spacing);
m_vtkVolumeImageData->GetOrigin(origin);
double center[3];
center[0] = origin[0] + spacing[0] * 0.5 * (extent[0] + extent[1]);
center[1] = origin[1] + spacing[1] * 0.5 * (extent[2] + extent[3]);
center[2] = origin[2] + spacing[2] * 0.5 * (extent[4] + extent[5]);
// Matrices for axial, coronal, sagittal, oblique view orientations
static double axialElements[16] = {
1, 0, 0, 0,
0, 1, 0, 0,
0, 0, 1, 0,
0, 0, 0, 1 };
static double coronalElements[16] = {
1, 0, 0, 0,
0, 0, 1, 0,
0,-1, 0, 0,
0, 0, 0, 1 };
static double sagittalElements[16] = {
0, 0,-1, 0,
1, 0, 0, 0,
0,-1, 0, 0,
0, 0, 0, 1 };
static double obliqueElements[16] = {
1, 0, 0, 0,
0, 0.866025, -0.5, 0,
0, 0.5, 0.866025, 0,
0, 0, 0, 1 };
// Set the slice orientation
vtkSmartPointer<vtkMatrix4x4> resliceAxes = vtkSmartPointer<vtkMatrix4x4>::New();
resliceAxes->DeepCopy(axialElements);
// Set the point through which to slice
resliceAxes->SetElement(0, 3, 0);
resliceAxes->SetElement(1, 3, 0);
resliceAxes->SetElement(2, 3, 70);
// Extract a slice in the desired orientation
vtkSmartPointer<vtkImageReslice> reslice = vtkSmartPointer<vtkImageReslice>::New();
reslice->SetInputData(m_vtkVolumeImageData);
reslice->SetOutputDimensionality(2);
reslice->SetResliceAxes(resliceAxes);
reslice->SetInterpolationModeToLinear();
reslice->Update();
// Create a lookup table
table = vtkSmartPointer<vtkLookupTable>::New();
// image intensity range
// my image is USHORT TYPE, so I Set 0 ~ 65535
table->SetRange(0, 65535);
table->SetNumberOfTableValues(65536);
double red = 0, green = 0, blue = 0;
for (int i = 0; i < table->GetNumberOfTableValues(); i++)
{
// 0~19999 value have Black ~ Red
if (i < 20000)
red += 1.0 / 20000;
// 20000~39999 value have Red ~ Yellow
else if (i < 40000)
green += 1.0 / 20000;
// 40000~59999 value have Yellow ~ white
// and 60000~ have white
else if (i < 60000)
blue += 1.0 / 20000;
table->SetTableValue(i, red, green, blue, 1);
}
table->Build();
//// Map the image through the lookup table
vtkSmartPointer<vtkImageMapToColors> color = vtkSmartPointer<vtkImageMapToColors>::New();
color->SetLookupTable(table);
color->SetInputData(reslice->GetOutput());
color->Update();
vtkSmartPointer<vtkDataSetMapper> mapper = vtkSmartPointer<vtkDataSetMapper>::New();
mapper->SetInputData(color->GetOutput());
//Setting ResliceImagePlaneActor
//Location, Opacity, Connect Mapper and Actor
double position[3] = { 0, 0, 70 };
vtkSmartPointer<vtkActor> nomal_actor = vtkSmartPointer<vtkActor>::New();
nomal_actor->SetMapper(mapper);
nomal_actor->GetProperty()->SetOpacity(0.7);
nomal_actor->SetPosition(position);
m_vtkRenderer->AddActor(nomal_actor);
//Setting ReslicePlaneOutLineActor
vtkSmartPointer<vtkOutlineFilter> sliceOutlineFilter = vtkSmartPointer<vtkOutlineFilter>::New();
sliceOutlineFilter->SetInputData(color->GetOutput());
sliceOutlineFilter->Update();
vtkSmartPointer<vtkPolyDataMapper> sliceOutlineMapper = vtkSmartPointer<vtkPolyDataMapper>::New();
sliceOutlineMapper->SetInputData(sliceOutlineFilter->GetOutput());
vtkSmartPointer<vtkActor> vtksliceOutlineActor = vtkSmartPointer<vtkActor>::New();
vtksliceOutlineActor->SetMapper(sliceOutlineMapper);
vtksliceOutlineActor->GetProperty()->SetColor(1, 0, 0);
vtksliceOutlineActor->SetPosition(position);
m_vtkRenderer->AddActor(vtksliceOutlineActor);
//Setting ScalarBarActor
//It is fine to Skip
vtkSmartPointer<vtkScalarBarActor> scalarBar = vtkSmartPointer<vtkScalarBarActor>::New();
scalarBar->SetLookupTable(table);
scalarBar->SetTitle("value");
scalarBar->SetNumberOfLabels(10);
scalarBar->SetLabelFormat("%10.2f");
scalarBar->SetHeight(.2);
scalarBar->SetWidth(.2);
m_vtkRenderer->AddActor2D(scalarBar);
}
这是我的完整代码
检查一下,如果您不理解某些部分,请Post评论:)
给定一个 DICOM 切片,位置直接进入 SetResliceAxesOrigin
,而纹理基元被移动一点以使像素单元的 VTK 中心与 OpenGL 边界相对应。然后沿 X 轴和 Y 轴的图像尺寸和间距用于计算图元的范围(跟踪 DICOM 标签的打印)。目前可以省略路线,因为 vtkImageData
没有路线 (https://discourse.vtk.org/t/proposal-to-add-orientation-to-vtkimagedata-feedback-wanted/120/2)。但它们很可能(希望)在阅读过程中被烘焙到图像数据中。
您可以通过与体积或表面模型进行比较来验证渲染。
对于处理医学图像数据,我推荐ITK (https://itk.org/)。如果仍然不使用 ITK,只需读取图像并在代码中将输出称为 image
。
打印样例:
0020|0037 1\-2.58172e-10.63164e-11\-1.26711e-11\-0.187259\-0.982311
Warning, direction will be lost
0018|0050 1
0028|0030 1
0020|0032 -252.69\-118.2735.807
#include "itkImage.h"
#include "itkImageFileReader.h"
#include "vtkInteractorStyleTrackballCamera.h"
#include "vtkImplicitPlaneWidget2.h"
#include "vtkImplicitPlaneRepresentation.h"
#include "vtkCommand.h"
#include "vtkPlane.h"
#include "vtkImageMapToColors.h"
#include "vtkSmartPointer.h"
#include "vtkPolyData.h"
#include "vtkPoints.h"
#include "vtkCellArray.h"
#include "vtkFloatArray.h"
#include "vtkImageReslice.h"
#include "vtkPiecewiseFunction.h"
#include "vtkColorTransferFunction.h"
#include "vtkGPUVolumeRayCastMapper.h"
#include "vtkVolumeProperty.h"
#include "vtkTexture.h"
#include "vtkLookupTable.h"
#include "vtkImageMarchingCubes.h"
vtkSmartPointer<vtkPolyData> MyCreateSimplePlane( const double * corners )
{
vtkSmartPointer<vtkPolyData> ret=vtkSmartPointer<vtkPolyData>::New();
vtkSmartPointer< vtkFloatArray > tcoords = vtkSmartPointer< vtkFloatArray >::New();
tcoords->SetNumberOfComponents( 2 );
tcoords->InsertNextTuple2( 0, 0 );
tcoords->InsertNextTuple2( 1, 0 );
tcoords->InsertNextTuple2( 1, 1 );
tcoords->InsertNextTuple2( 0, 1 );
ret->GetPointData()->SetTCoords( tcoords );
vtkSmartPointer<vtkFloatArray> floatArray = vtkSmartPointer<vtkFloatArray>::New();
floatArray->SetNumberOfComponents( 3 );
for ( int i=0;i<4;++i )
{
floatArray->InsertNextTuple3( corners[3*i], corners[3*i+1], corners[3*i+2] );
}
vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New();
points->SetData( floatArray );
ret->SetPoints( points );
vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
cells->InsertNextCell( 4 );
cells->InsertCellPoint( 0 );
cells->InsertCellPoint( 1 );
cells->InsertCellPoint( 2 );
cells->InsertCellPoint( 3 );
ret->SetPolys( cells );
return ret;
}
void mainactual( void )
{
bool doVolumeRender = true;
bool doSurfaceRender = false;
typedef itk::Image<short,3> Image;
typedef itk::VTKImageExport< Image > Export;
typedef itk::ImageFileReader<Image> Reader;
// read image using ITK (this is essentially an nrrd-converted Visible Human CT "Male")
// https://mri.radiology.uiowa.edu/VHDicom/
// prefer using ITK for medical data
Reader::Pointer reader = Reader::New();
reader->SetFileName("D:/tmp/vhp_male.nrrd");
reader->Update();
// if input image orientation is not (1,0,0[=11=],1,0), user should handle
// the respective transform by using, e.g., vtkVolume::SetOrientation
Image::DirectionType dir = reader->GetOutput()->GetDirection();
std::cout << "0020|0037" << " ";
std::cout << dir(0,0) << "\" << dir(1,0) << "\" << dir(2,0) << "\";
std::cout << dir(0,1) << "\" << dir(1,1) << "\" << dir(2,1);
std::cout << std::endl;
if ( !dir.GetVnlMatrix().is_identity( 0.0001 ) )
std::cout << "Warning, direction will be lost" << std::endl;
// import ITK image into VTK
vtkSmartPointer<vtkImageImport> imp = vtkSmartPointer<vtkImageImport>::New();
Export::Pointer exp = Export::New();
exp->SetInput( reader->GetOutput() );
imp->SetUpdateInformationCallback(exp->GetUpdateInformationCallback());
imp->SetPipelineModifiedCallback(exp->GetPipelineModifiedCallback());
imp->SetWholeExtentCallback(exp->GetWholeExtentCallback());
imp->SetSpacingCallback(exp->GetSpacingCallback());
imp->SetOriginCallback(exp->GetOriginCallback());
imp->SetScalarTypeCallback(exp->GetScalarTypeCallback());
imp->SetNumberOfComponentsCallback(exp->GetNumberOfComponentsCallback());
imp->SetPropagateUpdateExtentCallback(exp->GetPropagateUpdateExtentCallback());
imp->SetUpdateDataCallback(exp->GetUpdateDataCallback());
imp->SetDataExtentCallback(exp->GetDataExtentCallback());
imp->SetBufferPointerCallback(exp->GetBufferPointerCallback());
imp->SetCallbackUserData(exp->GetCallbackUserData());
imp->Update();
// this is our vtkImageData counterpart as it was read using VTK
vtkImageData * image = imp->GetOutput();
// create render window
vtkSmartPointer< vtkRenderer > ren = vtkSmartPointer< vtkRenderer >::New();
vtkSmartPointer< vtkRenderWindow > rw = vtkSmartPointer< vtkRenderWindow >::New();
rw->AddRenderer( ren );
rw->SetSize( 1024,1024 );
// create interactor used later
vtkSmartPointer< vtkRenderWindowInteractor > ia = vtkSmartPointer<vtkRenderWindowInteractor>::New();
ia->SetRenderWindow( rw );
ia->SetInteractorStyle( vtkSmartPointer< vtkInteractorStyleTrackballCamera >::New() );
ia->Initialize();
// define cutplane early on
vtkSmartPointer<vtkPlane > cutplane = vtkSmartPointer<vtkPlane>::New();
cutplane->SetNormal( 0,0,-1);
if ( doVolumeRender )
{
// set up some fancy volume rendering transfer function
vtkSmartPointer<vtkPiecewiseFunction> pw = vtkSmartPointer<vtkPiecewiseFunction>::New();
pw->AddPoint(-761.61130742049477, 0);
pw->AddPoint(-40.042402826855096, 0);
pw->AddPoint(353.54063604240287, 0.28817733990147787);
pw->AddPoint(1091.5088339222616, 0.69458128078817727);
pw->AddPoint(1763.8798586572439, 1);
vtkSmartPointer<vtkPiecewiseFunction> gf = vtkSmartPointer<vtkPiecewiseFunction>::New();
gf->AddPoint(-1024, 0);
gf->AddPoint(-1007.6007067137809, 1);
gf->AddPoint(-220.43462897526501, 0.78244274809160308);
gf->AddPoint(697.92579505300364, 0.9007633587786259);
gf->AddPoint(1157.1060070671379, 0.53435114503816794);
vtkSmartPointer<vtkColorTransferFunction> cf = vtkSmartPointer<vtkColorTransferFunction>::New();
cf->AddRGBPoint(-105.63957597173146, 0, 0, 0, 0.5, 0);
cf->AddRGBPoint(255.14487632508826, 0.93333299999999997, 0, 0, 0.5, 0);
cf->AddRGBPoint(353.54063604240287, 1, 0.90588199999999997, 0.66666700000000001, 0.5, 0);
cf->AddRGBPoint(632.32862190812716, 1, 0.66666666666666663, 0, 0.5, 0);
cf->AddRGBPoint(779.92226148409895, 1, 1, 1, 0.5, 0);
// and make GPUVolumeRayCastMapper to render them
typedef vtkGPUVolumeRayCastMapper Mapper;
vtkSmartPointer< Mapper > mapper = vtkSmartPointer< Mapper >::New();
mapper->SetInputData( image );
vtkSmartPointer< vtkVolumeProperty > prop = vtkSmartPointer< vtkVolumeProperty >::New();
prop->SetColor( cf );
prop->SetScalarOpacity( pw );
prop->SetGradientOpacity( gf );
prop->SetInterpolationTypeToLinear();
prop->SetDiffuse( 0 );
prop->SetSpecular( 0.0 );
prop->SetAmbient( 1 );
mapper->SetUseDepthPass( 1 );
mapper->SetUseJittering( 1 );
mapper->SetAutoAdjustSampleDistances( 0 );
vtkSmartPointer< vtkVolume > volume = vtkSmartPointer< vtkVolume >::New();
volume->SetMapper( mapper );
volume->SetProperty( prop );
// clip the volume
mapper->AddClippingPlane( cutplane );
ren->AddVolume( volume );
}
if ( doSurfaceRender )
{
// do marching cubes polygons
vtkSmartPointer<vtkImageMarchingCubes> mc = vtkSmartPointer<vtkImageMarchingCubes>::New();
mc->SetComputeGradients( 0 );
mc->SetComputeNormals( 0 );
mc->SetComputeScalars( 0 );
mc->SetInputData( image );
mc->SetNumberOfContours( 1 );
mc->SetValue( 0, 100.0 );
mc->Update();
vtkSmartPointer< vtkPolyDataMapper > surfmapper = vtkSmartPointer< vtkPolyDataMapper >::New();
vtkSmartPointer< vtkActor > surf = vtkSmartPointer< vtkActor >::New();
surf->SetMapper( surfmapper );
surfmapper->SetInputData( mc->GetOutput() );
surf->GetProperty()->SetAmbient(0);
surf->GetProperty()->SetDiffuse(1);
surf->GetProperty()->SetSpecular(0);
surf->GetProperty()->SetColor( 0.5, 0.9, 0.1 );
surfmapper->AddClippingPlane( cutplane );
ren->AddActor( surf );
}
// do the image slice plane
vtkSmartPointer< vtkImageReslice > slicer = vtkSmartPointer< vtkImageReslice >::New();
slicer->SetInputData( image );
vtkSmartPointer<vtkMatrix4x4> identity = vtkSmartPointer<vtkMatrix4x4>::New() ;
//identity->Identity(); // not needed as the orientation is identity anyways
slicer->SetResliceAxes( identity );
slicer->SetOutputDimensionality( 2 );
slicer->SetInterpolationModeToLinear();
vtkSmartPointer< vtkTexture > texture = vtkSmartPointer< vtkTexture > ::New();
vtkSmartPointer< vtkLookupTable > table = vtkSmartPointer< vtkLookupTable >::New();
table->SetNumberOfColors( 256 );
table->SetHueRange( 0.0, 0.0 );
table->SetSaturationRange( 0, 0 );
table->SetValueRange( 0.0, 1.0 );
table->SetAlphaRange( 1.0, 1.0 );
table->SetUseBelowRangeColor( 1 );
table->SetBelowRangeColor(0,0,0,0);
table->SetRange( -200, 200 );
table->Build();
texture->SetInputConnection( slicer->GetOutputPort() );
texture->SetLookupTable( table );
texture->SetColorModeToMapScalars();
vtkSmartPointer< vtkPolyDataMapper > polymapper = vtkSmartPointer< vtkPolyDataMapper >::New();
vtkSmartPointer< vtkActor > plane = vtkSmartPointer< vtkActor >::New();
plane->SetMapper( polymapper );
plane->GetProperty()->SetTexture(0, texture );
plane->GetProperty()->SetAmbient(1);
plane->GetProperty()->SetDiffuse(0);
plane->GetProperty()->SetSpecular(0);
plane->GetProperty()->SetAmbientColor(1,1,1);
plane->GetProperty()->SetEdgeColor(1,0,0);
plane->GetProperty()->SetEdgeVisibility(1);
ren->AddActor( plane );
int extent[6];
double origin[3];
double spacing[3];
image->GetOrigin( origin );
image->GetSpacing( spacing );
image->GetExtent( extent );
{
std::stringstream s;
std::cout << "0018|0050" << " ";
std::cout << spacing[2];
std::cout << std::endl;
}
{
std::cout << "0028|0030" << " ";
std::cout << spacing[0] << "\" << spacing[1];
std::cout << std::endl;
}
ren->ResetCamera(); // before looping, hope there is either surface or volume
int z = (extent[5]-extent[4])/2;
// for ( ... )
{
// DICOM pixel corner
double corner[3]={ origin[0], origin[1], origin[2] + (double)z * spacing[2] };
std::cout << "0020|0032" << " ";
std::cout << corner[0] << "\" << corner[1] << "\" << corner[2];
std::cout << std::endl;
slicer->SetResliceAxesOrigin( corner );
double corners[12];
// remove the half spacing to make DICOM / ITK / VTK origin to OpenGL origin
corners[0]=corner[0]-0.5*spacing[0];
corners[1]=corner[1]-0.5*spacing[1];
corners[2]=corner[2]; // the slicing direction
// +1 are to convert from DICOM pixel coordinates to OpenGL bounds
const double vx[3]={ (double)(extent[1]-extent[0]+1)*spacing[0], 0.0, 0.0 };
const double vy[3]={ 0.0, (double)(extent[3]-extent[2]+1)*spacing[1], 0.0} ;
for ( unsigned int u=0;u<3;++u )
{
corners[ 3+ u ] = corners[u] + vx[u];
corners[ 6+ u ] = corners[u] + vx[u] + vy[u];
corners[ 9+ u ] = corners[u] + vy[u];
}
cutplane->SetOrigin( corner );
vtkSmartPointer< vtkPolyData > polys = MyCreateSimplePlane( corners );
polymapper->SetInputData( polys );
ia->Render();
}
ia->Start();
}
* 更新 *
注意:方向必须单独处理,例如,通过 SetUserTransform
用于平面和体积,或图像数据的重新采样。下面显示的两个图像沿最外层图像维度进行类似切片,而不是沿世界空间 z 轴:
我正在做一个 VTK 程序,当我输入 2D DICOM 患者图像位置 (请参考给定的图像以便更好地理解),我需要在 3D 表面渲染输出中获取特定切片。
对于体积渲染的 3D 图像,这可以通过使用这些函数来实现,即 vtkImageData、vtkImageMapToColors、vtkImage Actor。
为了便于理解,我展示了一张示例图片
考虑一下我的 3d 输出将如下图所示
和当我在文本框中输入图像位置(患者)并单击确定按钮时,对应的切片应显示在 3d 图像中,如下图所示
如果我的问题无法理解,请告知 我被困在这里我什至不知道我的代码是 right.Here 是不是我的代码
ExtractVOI->SetInput(reader1->GetOutput());//VOI extractor
ExtractVOI->SetVOI(1,0,0,0,1,0);//i have given the Image Orientation(patient) as the SetVOI value
////====CREATE LookUpTable
tableax1->SetTableRange(0.0, 4096.0);
tableax1->SetValueRange(0.0, 1.0);
tableax1->SetSaturationRange(0.0, 0.0);
tableax1->SetRampToSCurve();
tableax1->SetAlphaRange(0.0, 0.08);
tableax1->Build();
//====CREATE TEXTURE
planesourceax1->SetXResolution(1);
planesourceax1->SetYResolution(1);
planesourceax1->SetOrigin(0,0,0);
planesourceax1->SetPoint1(xg , yg,zg);//i have given the value of Image Position(patient) that is taken from a textbox ,as the points
planesourceax1->Update();
vtkSmartPointer<vtkPolyDataMapper> mapax1 = vtkSmartPointer<vtkPolyDataMapper>::New();
mapax1->SetInputConnection(planesourceax1->GetOutputPort());
mapax1->UpdateWholeExtent();
textureax1->SetInputConnection(ExtractVOI->GetOutputPort());
textureax1->InterpolateOn();
textureax1->SetLookupTable(tableax1);
textureax1->UpdateWholeExtent();
//===PASS TO ACTOR
actorax1->SetMapper(mapax1);
actorax1->GetMapper()->SetResolveCoincidentTopologyToPolygonOffset();
actorax1->GetMapper()->SetResolveCoincidentTopologyPolygonOffsetParameters(0.1, -1.0);
actorax1->SetTexture(textureax1);
renderer->AddActor(actorax1);
renderWindow->Render();
但我没有得到输出
我也试过了:
static double axialElements[16] = {
1, 0, 0, 0,
0, 1, 0, 0,
0, 0, 1, 0,
0, 0, 0, 1 };
resliceax1->SetInputConnection(reader->GetOutputPort());
resliceax1->SetOutputDimensionality(2);
vtkSmartPointer<vtkMatrix4x4> reslicematrixax1 = vtkSmartPointer<vtkMatrix4x4>::New();
reslicematrixax1->DeepCopy(axialElements);
resliceax1->SetResliceAxes(reslicematrixax1);
resliceax1->SetResliceAxesOrigin(0.0, 0.0, 0.0);
resliceax1->Update();
extractaxpos1->RemoveAllInputs();
extractaxpos1->SetInputConnection(resliceax1->GetOutputPort());
////====CREATE LUT
tableax1->SetTableRange(0.0, 4096.0);
tableax1->SetValueRange(0.0, 1.0);
tableax1->SetSaturationRange(0.0, 0.0);
tableax1->SetRampToSCurve();
tableax1->SetAlphaRange(0.0, 0.08);
tableax1->Build();
//====CREATE TEXTURE
planesourceax1->SetXResolution(1);
planesourceax1->SetYResolution(1);
planesourceax1->SetOrigin(0,0,0);
planesourceax1->SetPoint1((xval/20 + xval/32),(yval/20 + yval/32),(zval/20 + zval/32));//this is where i put the values ad divided by its tag id(0020,0032)
//planesourceax1->SetPoint2(fBounds[0] , fBounds[3], fBounds[4]);
planesourceax1->Update();
vtkSmartPointer<vtkPolyDataMapper> mapax1 = vtkSmartPointer<vtkPolyDataMapper>::New();
mapax1->SetInputConnection(planesourceax1->GetOutputPort());
mapax1->UpdateWholeExtent();
textureax1->SetInputConnection(extractaxpos1->GetOutputPort());
textureax1->InterpolateOn();
textureax1->SetLookupTable(tableax1);
textureax1->UpdateWholeExtent();
//===PASS TO ACTOR
actorax1->SetMapper(mapax1);
actorax1->GetMapper()->SetResolveCoincidentTopologyToPolygonOffset();
actorax1->GetMapper()->SetResolveCoincidentTopologyPolygonOffsetParameters(0.1, -1.0);
actorax1->SetTexture(textureax1);
resliceax1->SetResliceAxesOrigin(0.0, 0.0, 0.0);
actorax1->SetPosition((xval/20 + xval/32),(yval/20 + yval/32),(zval/20 + zval/32));//I made the same changes here also
planesourceax1->SetOrigin(fBoundsUpdated[0], fBoundsUpdated[2], pDoc->fBounds[4]);
planesourceax1->SetPoint1(fBoundsUpdated[1] , fBoundsUpdated[2], pDoc->fBounds[4]);
planesourceax1->SetPoint2(fBoundsUpdated[0] , fBoundsUpdated[3], pDoc->fBounds[4]);
planesourceax1->Update();
但不是切割切片的位置is.It是切割不同的位置
如果你想知道什么元素决定切割方向
https://github.com/Kitware/VTK/blob/master/Examples/ImageProcessing/Cxx/ImageSlicing.cxx
请看这个link
static double axialElements[16] = {
1, 0, 0, 0,
0, 1, 0, 0,
0, 0, 1, 0,
0, 0, 0, 1 };
...
vtkSmartPointer<vtkMatrix4x4> resliceAxes = vtkSmartPointer<vtkMatrix4x4>::New();
resliceAxes->DeepCopy(axialElements);
// Set the point through which to slice
resliceAxes->SetElement(0, 3, 0);
resliceAxes->SetElement(1, 3, 0);
resliceAxes->SetElement(2, 3, 70);
矩阵 4x4 正在设置 重新切片方向 和 SetElement 决定 Reslice origin
此代码设置选项为Reslice XYimage start at Z = 70
vtkSmartPointer<vtkImageReslice> reslice = vtkSmartPointer<vtkImageReslice>::New();
reslice->SetInputData(m_vtkVolumeImageData);
reslice->SetOutputDimensionality(2);
reslice->SetResliceAxes(resliceAxes);
reslice->SetInterpolationModeToLinear();
reslice->Update();
此代码生成 Reslice 数据。 结果维度为 2D(平面) 如果你设置这个顺序, 将 Reslice 数据连接到 vtkImageMapToColors 以绘制 Reslice 图像。
最后,连接 Mapper 和 Actor 以使用 Volume 进行展示。 您需要设置 Actor 位置,因为 Reslice Data 可能没有位置信息
如果平面位置没有改变,使用 vtkDataSetMapper 和 vtkActor,而不是 vtkImageActor
或者你只是使用Widget, 我推荐 vtkImagePlaneWidget。
操作简单,功能强大。
希望对您有所帮助。
如果需要完整代码,请告诉我
void CMFCApplication1Dlg::setSliceImage()
{
int extent[6];
double spacing[3];
double origin[3];
//m_vtkVolumeImageData is vtkImageData(Save Slice Images)
//Same as dicomReader->GetOutput();
m_vtkVolumeImageData->GetExtent(extent);
m_vtkVolumeImageData->GetSpacing(spacing);
m_vtkVolumeImageData->GetOrigin(origin);
double center[3];
center[0] = origin[0] + spacing[0] * 0.5 * (extent[0] + extent[1]);
center[1] = origin[1] + spacing[1] * 0.5 * (extent[2] + extent[3]);
center[2] = origin[2] + spacing[2] * 0.5 * (extent[4] + extent[5]);
// Matrices for axial, coronal, sagittal, oblique view orientations
static double axialElements[16] = {
1, 0, 0, 0,
0, 1, 0, 0,
0, 0, 1, 0,
0, 0, 0, 1 };
static double coronalElements[16] = {
1, 0, 0, 0,
0, 0, 1, 0,
0,-1, 0, 0,
0, 0, 0, 1 };
static double sagittalElements[16] = {
0, 0,-1, 0,
1, 0, 0, 0,
0,-1, 0, 0,
0, 0, 0, 1 };
static double obliqueElements[16] = {
1, 0, 0, 0,
0, 0.866025, -0.5, 0,
0, 0.5, 0.866025, 0,
0, 0, 0, 1 };
// Set the slice orientation
vtkSmartPointer<vtkMatrix4x4> resliceAxes = vtkSmartPointer<vtkMatrix4x4>::New();
resliceAxes->DeepCopy(axialElements);
// Set the point through which to slice
resliceAxes->SetElement(0, 3, 0);
resliceAxes->SetElement(1, 3, 0);
resliceAxes->SetElement(2, 3, 70);
// Extract a slice in the desired orientation
vtkSmartPointer<vtkImageReslice> reslice = vtkSmartPointer<vtkImageReslice>::New();
reslice->SetInputData(m_vtkVolumeImageData);
reslice->SetOutputDimensionality(2);
reslice->SetResliceAxes(resliceAxes);
reslice->SetInterpolationModeToLinear();
reslice->Update();
// Create a lookup table
table = vtkSmartPointer<vtkLookupTable>::New();
// image intensity range
// my image is USHORT TYPE, so I Set 0 ~ 65535
table->SetRange(0, 65535);
table->SetNumberOfTableValues(65536);
double red = 0, green = 0, blue = 0;
for (int i = 0; i < table->GetNumberOfTableValues(); i++)
{
// 0~19999 value have Black ~ Red
if (i < 20000)
red += 1.0 / 20000;
// 20000~39999 value have Red ~ Yellow
else if (i < 40000)
green += 1.0 / 20000;
// 40000~59999 value have Yellow ~ white
// and 60000~ have white
else if (i < 60000)
blue += 1.0 / 20000;
table->SetTableValue(i, red, green, blue, 1);
}
table->Build();
//// Map the image through the lookup table
vtkSmartPointer<vtkImageMapToColors> color = vtkSmartPointer<vtkImageMapToColors>::New();
color->SetLookupTable(table);
color->SetInputData(reslice->GetOutput());
color->Update();
vtkSmartPointer<vtkDataSetMapper> mapper = vtkSmartPointer<vtkDataSetMapper>::New();
mapper->SetInputData(color->GetOutput());
//Setting ResliceImagePlaneActor
//Location, Opacity, Connect Mapper and Actor
double position[3] = { 0, 0, 70 };
vtkSmartPointer<vtkActor> nomal_actor = vtkSmartPointer<vtkActor>::New();
nomal_actor->SetMapper(mapper);
nomal_actor->GetProperty()->SetOpacity(0.7);
nomal_actor->SetPosition(position);
m_vtkRenderer->AddActor(nomal_actor);
//Setting ReslicePlaneOutLineActor
vtkSmartPointer<vtkOutlineFilter> sliceOutlineFilter = vtkSmartPointer<vtkOutlineFilter>::New();
sliceOutlineFilter->SetInputData(color->GetOutput());
sliceOutlineFilter->Update();
vtkSmartPointer<vtkPolyDataMapper> sliceOutlineMapper = vtkSmartPointer<vtkPolyDataMapper>::New();
sliceOutlineMapper->SetInputData(sliceOutlineFilter->GetOutput());
vtkSmartPointer<vtkActor> vtksliceOutlineActor = vtkSmartPointer<vtkActor>::New();
vtksliceOutlineActor->SetMapper(sliceOutlineMapper);
vtksliceOutlineActor->GetProperty()->SetColor(1, 0, 0);
vtksliceOutlineActor->SetPosition(position);
m_vtkRenderer->AddActor(vtksliceOutlineActor);
//Setting ScalarBarActor
//It is fine to Skip
vtkSmartPointer<vtkScalarBarActor> scalarBar = vtkSmartPointer<vtkScalarBarActor>::New();
scalarBar->SetLookupTable(table);
scalarBar->SetTitle("value");
scalarBar->SetNumberOfLabels(10);
scalarBar->SetLabelFormat("%10.2f");
scalarBar->SetHeight(.2);
scalarBar->SetWidth(.2);
m_vtkRenderer->AddActor2D(scalarBar);
}
这是我的完整代码
检查一下,如果您不理解某些部分,请Post评论:)
给定一个 DICOM 切片,位置直接进入 SetResliceAxesOrigin
,而纹理基元被移动一点以使像素单元的 VTK 中心与 OpenGL 边界相对应。然后沿 X 轴和 Y 轴的图像尺寸和间距用于计算图元的范围(跟踪 DICOM 标签的打印)。目前可以省略路线,因为 vtkImageData
没有路线 (https://discourse.vtk.org/t/proposal-to-add-orientation-to-vtkimagedata-feedback-wanted/120/2)。但它们很可能(希望)在阅读过程中被烘焙到图像数据中。
您可以通过与体积或表面模型进行比较来验证渲染。
对于处理医学图像数据,我推荐ITK (https://itk.org/)。如果仍然不使用 ITK,只需读取图像并在代码中将输出称为 image
。
打印样例:
0020|0037 1\-2.58172e-10.63164e-11\-1.26711e-11\-0.187259\-0.982311
Warning, direction will be lost
0018|0050 1
0028|0030 1
0020|0032 -252.69\-118.2735.807
#include "itkImage.h"
#include "itkImageFileReader.h"
#include "vtkInteractorStyleTrackballCamera.h"
#include "vtkImplicitPlaneWidget2.h"
#include "vtkImplicitPlaneRepresentation.h"
#include "vtkCommand.h"
#include "vtkPlane.h"
#include "vtkImageMapToColors.h"
#include "vtkSmartPointer.h"
#include "vtkPolyData.h"
#include "vtkPoints.h"
#include "vtkCellArray.h"
#include "vtkFloatArray.h"
#include "vtkImageReslice.h"
#include "vtkPiecewiseFunction.h"
#include "vtkColorTransferFunction.h"
#include "vtkGPUVolumeRayCastMapper.h"
#include "vtkVolumeProperty.h"
#include "vtkTexture.h"
#include "vtkLookupTable.h"
#include "vtkImageMarchingCubes.h"
vtkSmartPointer<vtkPolyData> MyCreateSimplePlane( const double * corners )
{
vtkSmartPointer<vtkPolyData> ret=vtkSmartPointer<vtkPolyData>::New();
vtkSmartPointer< vtkFloatArray > tcoords = vtkSmartPointer< vtkFloatArray >::New();
tcoords->SetNumberOfComponents( 2 );
tcoords->InsertNextTuple2( 0, 0 );
tcoords->InsertNextTuple2( 1, 0 );
tcoords->InsertNextTuple2( 1, 1 );
tcoords->InsertNextTuple2( 0, 1 );
ret->GetPointData()->SetTCoords( tcoords );
vtkSmartPointer<vtkFloatArray> floatArray = vtkSmartPointer<vtkFloatArray>::New();
floatArray->SetNumberOfComponents( 3 );
for ( int i=0;i<4;++i )
{
floatArray->InsertNextTuple3( corners[3*i], corners[3*i+1], corners[3*i+2] );
}
vtkSmartPointer<vtkPoints> points = vtkSmartPointer<vtkPoints>::New();
points->SetData( floatArray );
ret->SetPoints( points );
vtkSmartPointer<vtkCellArray> cells = vtkSmartPointer<vtkCellArray>::New();
cells->InsertNextCell( 4 );
cells->InsertCellPoint( 0 );
cells->InsertCellPoint( 1 );
cells->InsertCellPoint( 2 );
cells->InsertCellPoint( 3 );
ret->SetPolys( cells );
return ret;
}
void mainactual( void )
{
bool doVolumeRender = true;
bool doSurfaceRender = false;
typedef itk::Image<short,3> Image;
typedef itk::VTKImageExport< Image > Export;
typedef itk::ImageFileReader<Image> Reader;
// read image using ITK (this is essentially an nrrd-converted Visible Human CT "Male")
// https://mri.radiology.uiowa.edu/VHDicom/
// prefer using ITK for medical data
Reader::Pointer reader = Reader::New();
reader->SetFileName("D:/tmp/vhp_male.nrrd");
reader->Update();
// if input image orientation is not (1,0,0[=11=],1,0), user should handle
// the respective transform by using, e.g., vtkVolume::SetOrientation
Image::DirectionType dir = reader->GetOutput()->GetDirection();
std::cout << "0020|0037" << " ";
std::cout << dir(0,0) << "\" << dir(1,0) << "\" << dir(2,0) << "\";
std::cout << dir(0,1) << "\" << dir(1,1) << "\" << dir(2,1);
std::cout << std::endl;
if ( !dir.GetVnlMatrix().is_identity( 0.0001 ) )
std::cout << "Warning, direction will be lost" << std::endl;
// import ITK image into VTK
vtkSmartPointer<vtkImageImport> imp = vtkSmartPointer<vtkImageImport>::New();
Export::Pointer exp = Export::New();
exp->SetInput( reader->GetOutput() );
imp->SetUpdateInformationCallback(exp->GetUpdateInformationCallback());
imp->SetPipelineModifiedCallback(exp->GetPipelineModifiedCallback());
imp->SetWholeExtentCallback(exp->GetWholeExtentCallback());
imp->SetSpacingCallback(exp->GetSpacingCallback());
imp->SetOriginCallback(exp->GetOriginCallback());
imp->SetScalarTypeCallback(exp->GetScalarTypeCallback());
imp->SetNumberOfComponentsCallback(exp->GetNumberOfComponentsCallback());
imp->SetPropagateUpdateExtentCallback(exp->GetPropagateUpdateExtentCallback());
imp->SetUpdateDataCallback(exp->GetUpdateDataCallback());
imp->SetDataExtentCallback(exp->GetDataExtentCallback());
imp->SetBufferPointerCallback(exp->GetBufferPointerCallback());
imp->SetCallbackUserData(exp->GetCallbackUserData());
imp->Update();
// this is our vtkImageData counterpart as it was read using VTK
vtkImageData * image = imp->GetOutput();
// create render window
vtkSmartPointer< vtkRenderer > ren = vtkSmartPointer< vtkRenderer >::New();
vtkSmartPointer< vtkRenderWindow > rw = vtkSmartPointer< vtkRenderWindow >::New();
rw->AddRenderer( ren );
rw->SetSize( 1024,1024 );
// create interactor used later
vtkSmartPointer< vtkRenderWindowInteractor > ia = vtkSmartPointer<vtkRenderWindowInteractor>::New();
ia->SetRenderWindow( rw );
ia->SetInteractorStyle( vtkSmartPointer< vtkInteractorStyleTrackballCamera >::New() );
ia->Initialize();
// define cutplane early on
vtkSmartPointer<vtkPlane > cutplane = vtkSmartPointer<vtkPlane>::New();
cutplane->SetNormal( 0,0,-1);
if ( doVolumeRender )
{
// set up some fancy volume rendering transfer function
vtkSmartPointer<vtkPiecewiseFunction> pw = vtkSmartPointer<vtkPiecewiseFunction>::New();
pw->AddPoint(-761.61130742049477, 0);
pw->AddPoint(-40.042402826855096, 0);
pw->AddPoint(353.54063604240287, 0.28817733990147787);
pw->AddPoint(1091.5088339222616, 0.69458128078817727);
pw->AddPoint(1763.8798586572439, 1);
vtkSmartPointer<vtkPiecewiseFunction> gf = vtkSmartPointer<vtkPiecewiseFunction>::New();
gf->AddPoint(-1024, 0);
gf->AddPoint(-1007.6007067137809, 1);
gf->AddPoint(-220.43462897526501, 0.78244274809160308);
gf->AddPoint(697.92579505300364, 0.9007633587786259);
gf->AddPoint(1157.1060070671379, 0.53435114503816794);
vtkSmartPointer<vtkColorTransferFunction> cf = vtkSmartPointer<vtkColorTransferFunction>::New();
cf->AddRGBPoint(-105.63957597173146, 0, 0, 0, 0.5, 0);
cf->AddRGBPoint(255.14487632508826, 0.93333299999999997, 0, 0, 0.5, 0);
cf->AddRGBPoint(353.54063604240287, 1, 0.90588199999999997, 0.66666700000000001, 0.5, 0);
cf->AddRGBPoint(632.32862190812716, 1, 0.66666666666666663, 0, 0.5, 0);
cf->AddRGBPoint(779.92226148409895, 1, 1, 1, 0.5, 0);
// and make GPUVolumeRayCastMapper to render them
typedef vtkGPUVolumeRayCastMapper Mapper;
vtkSmartPointer< Mapper > mapper = vtkSmartPointer< Mapper >::New();
mapper->SetInputData( image );
vtkSmartPointer< vtkVolumeProperty > prop = vtkSmartPointer< vtkVolumeProperty >::New();
prop->SetColor( cf );
prop->SetScalarOpacity( pw );
prop->SetGradientOpacity( gf );
prop->SetInterpolationTypeToLinear();
prop->SetDiffuse( 0 );
prop->SetSpecular( 0.0 );
prop->SetAmbient( 1 );
mapper->SetUseDepthPass( 1 );
mapper->SetUseJittering( 1 );
mapper->SetAutoAdjustSampleDistances( 0 );
vtkSmartPointer< vtkVolume > volume = vtkSmartPointer< vtkVolume >::New();
volume->SetMapper( mapper );
volume->SetProperty( prop );
// clip the volume
mapper->AddClippingPlane( cutplane );
ren->AddVolume( volume );
}
if ( doSurfaceRender )
{
// do marching cubes polygons
vtkSmartPointer<vtkImageMarchingCubes> mc = vtkSmartPointer<vtkImageMarchingCubes>::New();
mc->SetComputeGradients( 0 );
mc->SetComputeNormals( 0 );
mc->SetComputeScalars( 0 );
mc->SetInputData( image );
mc->SetNumberOfContours( 1 );
mc->SetValue( 0, 100.0 );
mc->Update();
vtkSmartPointer< vtkPolyDataMapper > surfmapper = vtkSmartPointer< vtkPolyDataMapper >::New();
vtkSmartPointer< vtkActor > surf = vtkSmartPointer< vtkActor >::New();
surf->SetMapper( surfmapper );
surfmapper->SetInputData( mc->GetOutput() );
surf->GetProperty()->SetAmbient(0);
surf->GetProperty()->SetDiffuse(1);
surf->GetProperty()->SetSpecular(0);
surf->GetProperty()->SetColor( 0.5, 0.9, 0.1 );
surfmapper->AddClippingPlane( cutplane );
ren->AddActor( surf );
}
// do the image slice plane
vtkSmartPointer< vtkImageReslice > slicer = vtkSmartPointer< vtkImageReslice >::New();
slicer->SetInputData( image );
vtkSmartPointer<vtkMatrix4x4> identity = vtkSmartPointer<vtkMatrix4x4>::New() ;
//identity->Identity(); // not needed as the orientation is identity anyways
slicer->SetResliceAxes( identity );
slicer->SetOutputDimensionality( 2 );
slicer->SetInterpolationModeToLinear();
vtkSmartPointer< vtkTexture > texture = vtkSmartPointer< vtkTexture > ::New();
vtkSmartPointer< vtkLookupTable > table = vtkSmartPointer< vtkLookupTable >::New();
table->SetNumberOfColors( 256 );
table->SetHueRange( 0.0, 0.0 );
table->SetSaturationRange( 0, 0 );
table->SetValueRange( 0.0, 1.0 );
table->SetAlphaRange( 1.0, 1.0 );
table->SetUseBelowRangeColor( 1 );
table->SetBelowRangeColor(0,0,0,0);
table->SetRange( -200, 200 );
table->Build();
texture->SetInputConnection( slicer->GetOutputPort() );
texture->SetLookupTable( table );
texture->SetColorModeToMapScalars();
vtkSmartPointer< vtkPolyDataMapper > polymapper = vtkSmartPointer< vtkPolyDataMapper >::New();
vtkSmartPointer< vtkActor > plane = vtkSmartPointer< vtkActor >::New();
plane->SetMapper( polymapper );
plane->GetProperty()->SetTexture(0, texture );
plane->GetProperty()->SetAmbient(1);
plane->GetProperty()->SetDiffuse(0);
plane->GetProperty()->SetSpecular(0);
plane->GetProperty()->SetAmbientColor(1,1,1);
plane->GetProperty()->SetEdgeColor(1,0,0);
plane->GetProperty()->SetEdgeVisibility(1);
ren->AddActor( plane );
int extent[6];
double origin[3];
double spacing[3];
image->GetOrigin( origin );
image->GetSpacing( spacing );
image->GetExtent( extent );
{
std::stringstream s;
std::cout << "0018|0050" << " ";
std::cout << spacing[2];
std::cout << std::endl;
}
{
std::cout << "0028|0030" << " ";
std::cout << spacing[0] << "\" << spacing[1];
std::cout << std::endl;
}
ren->ResetCamera(); // before looping, hope there is either surface or volume
int z = (extent[5]-extent[4])/2;
// for ( ... )
{
// DICOM pixel corner
double corner[3]={ origin[0], origin[1], origin[2] + (double)z * spacing[2] };
std::cout << "0020|0032" << " ";
std::cout << corner[0] << "\" << corner[1] << "\" << corner[2];
std::cout << std::endl;
slicer->SetResliceAxesOrigin( corner );
double corners[12];
// remove the half spacing to make DICOM / ITK / VTK origin to OpenGL origin
corners[0]=corner[0]-0.5*spacing[0];
corners[1]=corner[1]-0.5*spacing[1];
corners[2]=corner[2]; // the slicing direction
// +1 are to convert from DICOM pixel coordinates to OpenGL bounds
const double vx[3]={ (double)(extent[1]-extent[0]+1)*spacing[0], 0.0, 0.0 };
const double vy[3]={ 0.0, (double)(extent[3]-extent[2]+1)*spacing[1], 0.0} ;
for ( unsigned int u=0;u<3;++u )
{
corners[ 3+ u ] = corners[u] + vx[u];
corners[ 6+ u ] = corners[u] + vx[u] + vy[u];
corners[ 9+ u ] = corners[u] + vy[u];
}
cutplane->SetOrigin( corner );
vtkSmartPointer< vtkPolyData > polys = MyCreateSimplePlane( corners );
polymapper->SetInputData( polys );
ia->Render();
}
ia->Start();
}
* 更新 *
注意:方向必须单独处理,例如,通过 SetUserTransform
用于平面和体积,或图像数据的重新采样。下面显示的两个图像沿最外层图像维度进行类似切片,而不是沿世界空间 z 轴: