对输入二维图像序列dicom进行二维图像分割后使用marching cube进行三维面绘制
在使用分割算法对二维图像序列进行分割后,将分割后的二维图像按照分割前一样的Marching cube算法就能重建得到需要的三维面数据了吗?我会实践之后贴上代码。
·
结合多篇论文,得到了常规流程:
在使用分割算法对二维图像序列进行分割后,将分割后的二维图像按照分割前一样的Marching cube算法就能重建得到需要的三维面数据。
#include <iostream>
#include <string>
#include <vector>
#include "itkImageSeriesReader.h"
#include "itkImageSeriesWriter.h"
#include "itkGDCMImageIO.h"
#include "itkGDCMSeriesFileNames.h"
#include "itkImageToVTKImageFilter.h"
#include "itkVTKImageToImageFilter.h"
#include "itkBinaryThresholdImageFilter.h"
#include "itkConnectedComponentImageFilter.h"
#include "itkRelabelComponentImageFilter.h"
#include "itkRescaleIntensityImageFilter.h"
#include "itkImageFileWriter.h"
#include "vtkSmartPointer.h"
#include "vtkPolyData.h"
#include "vtkPolyDataMapper.h"
#include "vtkActor.h"
#include "vtkRenderWindow.h"
#include "vtkRenderer.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkStructuredPoints.h"
#include "vtkStructuredPointsWriter.h"
typedef unsigned short PixelType;
const int Dimension = 3;
typedef itk::Image<PixelType, Dimension> ImageType;
typedef itk::ImageSeriesReader<ImageType> ReaderType;
typedef itk::ImageSeriesWriter<ImageType, itk::ImageFileWriter<ImageType>> WriterType;
typedef itk::GDCMImageIO ImageIOType;
typedef itk::GDCMSeriesFileNames NamesGeneratorType;
typedef itk::ImageToVTKImageFilter<ImageType> ImageToVTKFilterType;
typedef itk::VTKImageToImageFilter<ImageType> VTKToImageFilterType;
typedef itk::BinaryThresholdImageFilter<ImageType, ImageType> ThresholdFilterType;
typedef itk::ConnectedComponentImageFilter<ImageType, ImageType> ConnectedComponentFilterType;
typedef itk::RelabelComponentImageFilter<ImageType, ImageType> RelabelFilterType;
typedef itk::RescaleIntensityImageFilter<ImageType, ImageType> RescaleFilterType;
int main(int argc, char *argv[])
{
if (argc < 2) {
std::cerr << "Usage: " << argv[0] << " dicom_directory" << std::endl;
return EXIT_FAILURE;
}
// Step 1: Read DICOM series
std::string directory = argv[1];
ImageIOType::Pointer gdcmIO = ImageIOType::New();
NamesGeneratorType::Pointer namesGenerator = NamesGeneratorType::New();
namesGenerator->SetInputDirectory(directory);
const ReaderType::FileNamesContainer &fileNames = namesGenerator->GetInputFileNames();
ReaderType::Pointer reader = ReaderType::New();
reader->SetImageIO(gdcmIO);
reader->SetFileNames(fileNames);
try {
reader->Update();
}
catch (itk::ExceptionObject &ex) {
std::cerr << "Exception caught!" << std::endl;
std::cerr << ex << std::endl;
return EXIT_FAILURE;
}
// Step 2: Thresholding
ThresholdFilterType::Pointer thresholdFilter = ThresholdFilterType::New();
thresholdFilter->SetInput(reader->GetOutput());
thresholdFilter->SetLowerThreshold(200);
thresholdFilter->SetUpperThreshold(255);
thresholdFilter->SetInsideValue(1);
thresholdFilter->SetOutsideValue(0);
thresholdFilter->Update();
// Step 3: Connected component labeling
ConnectedComponentFilterType::Pointer connectedComponentFilter = ConnectedComponentFilterType::New();
connectedComponentFilter->SetInput(thresholdFilter->GetOutput());
connectedComponentFilter->Update();
// Step 4: Relabeling
RelabelFilterType::Pointer relabelFilter = RelabelFilterType::New();
relabelFilter->SetInput(connectedComponentFilter->GetOutput());
relabelFilter->Update();
// Step 5: Rescale intensity
RescaleFilterType::Pointer rescaleFilter = RescaleFilterType::New();
rescaleFilter->SetInput(relabelFilter->GetOutput());
rescaleFilter->SetOutputMinimum(0);
rescaleFilter->SetOutputMaximum(255);
rescaleFilter->Update();
// Step 6: Write output series
WriterType::Pointer writer = WriterType::New();
writer->SetInput(rescaleFilter->GetOutput());
writer->SetImageIO(gdcmIO);
writer->SetFileNames(namesGenerator->GetOutputFileNames());
writer->Update();
// Step 7: Convert to VTK image
ImageToVTKFilterType::Pointer imageToVTKFilter = ImageToVTKFilterType::New();
imageToVTKFilter->SetInput(rescaleFilter->GetOutput());
imageToVTKFilter->Update();
// Step 8: Convert to VTK polydata
vtkSmartPointer<vtkStructuredPoints> structuredPoints = vtkSmartPointer<vtkStructuredPoints>::New();
structuredPoints->SetDimensions(imageToVTKFilter->GetOutput()->GetDimensions());
structuredPoints->SetSpacing(imageToVTKFilter->GetOutput()->GetSpacing());
structuredPoints->SetOrigin(imageToVTKFilter->GetOutput()->GetOrigin());
structuredPoints->SetScalarTypeToUnsignedChar();
structuredPoints->GetPointData()->SetScalars(imageToVTKFilter->GetOutput()->GetPointData()->GetScalars());
vtkSmartPointer<vtkPolyData> polyData = vtkSmartPointer<vtkPolyData>::New();
vtkSmartPointer<vtkStructuredPointsWriter> writer2 = vtkSmartPointer<vtkStructuredPointsWriter>::New();
writer2->SetInputData(structuredPoints);
writer2->SetFileName("blood_vessels.vtk");
writer2->Write();
return EXIT_SUCCESS;
}
该程序主要分为以下几个步骤:
1.使用ITK库读取DICOM序列图像。
2.使用二值阈值滤波器将图像二值化。
3.使用连通成分标记滤波器标记和提取感兴趣区域。
4.使用重新标记滤波器对连通成分进行重新标记。
5.使用强度重缩放滤波器将像素值范围从[0,1]重缩放到[0,255]。
6.使用ITK库写入DICOM序列图像。
7.使用ITK-VTK转换器将ITK图像转换为VTK图像。
8.将VTK图像转换为VTK多边形数据,并使用VTK库将数据写入VTK文件。

魔乐社区(Modelers.cn) 是一个中立、公益的人工智能社区,提供人工智能工具、模型、数据的托管、展示与应用协同服务,为人工智能开发及爱好者搭建开放的学习交流平台。社区通过理事会方式运作,由全产业链共同建设、共同运营、共同享有,推动国产AI生态繁荣发展。
更多推荐
所有评论(0)