itk读取dicom,区域增长分割后vtk重建显示

int main(int argc, char *argv[])
{
  using InternalPixelType = float;
  constexpr unsigned int Dimension = 3;
  using InternalImageType = itk::Image<InternalPixelType, Dimension>;
  using OutputPixelType = float;
  using OutputImageType = itk::Image<OutputPixelType, Dimension>;
 
  using CastingFilterType =
      itk::CastImageFilter<InternalImageType, OutputImageType>;
  CastingFilterType::Pointer caster = CastingFilterType::New();
 
  using ReaderType = itk::ImageSeriesReader<InternalImageType>;
  using WriterType = itk::ImageFileWriter<OutputImageType>;
 
  using ImageIOType = itk::GDCMImageIO;
  using NamesGeneratorType = itk::GDCMSeriesFileNames;
  ReaderType::Pointer reader = ReaderType::New();
  WriterType::Pointer writer = WriterType::New();
 
  ImageIOType::Pointer dicomIO = ImageIOType::New();
  reader->SetImageIO(dicomIO);
  NamesGeneratorType::Pointer nameGenerator = NamesGeneratorType::New();
  nameGenerator->SetUseSeriesDetails(true);
  nameGenerator->SetDirectory("D:\\software\\30dicom\\222\\14");
  writer->SetFileName("D:\\software\\30dicom\\cttt.mha");
 
  using SeriesIdContainer = std::vector<std::string>;
  const SeriesIdContainer &seriesUID = nameGenerator->GetSeriesUIDs();
  auto seriesItr = seriesUID.begin();
  auto seriesEnd = seriesUID.end();
 
  using FileNamesContainer = std::vector<std::string>;
  FileNamesContainer fileNames;
  std::string seriesIdentifier;
  while (seriesItr != seriesEnd)
  {
 
    seriesIdentifier = seriesItr->c_str();
    fileNames = nameGenerator->GetFileNames(seriesIdentifier);
    ++seriesItr;
  }
 
  reader->SetFileNames(fileNames);
 
  try
  {
    reader->Update();
  }
  catch (itk::ExceptionObject &ex)
  {
    std::cout << ex << std::endl;
    return EXIT_FAILURE;
  }
  InternalImageType::SizeType imgSize = reader->GetOutput()->GetLargestPossibleRegion().GetSize();
  cout << "read done！Original size: " << imgSize << endl;
 
  //区域生长分割///
 
  using CurvatureFlowImageFilterType =
      itk::CurvatureFlowImageFilter<InternalImageType, InternalImageType>;
  CurvatureFlowImageFilterType::Pointer smoothing =
      CurvatureFlowImageFilterType::New();
 
  using ConnectedFilterType =
      itk::ConnectedThresholdImageFilter<InternalImageType, InternalImageType>;
  ConnectedFilterType::Pointer connectedThreshold =
      ConnectedFilterType::New();
 
  //smoothing->SetInput(reader->GetOutput());
  // smoothing->SetNumberOfIterations(5);
  // smoothing->SetTimeStep(0.125);
  connectedThreshold->SetInput(reader->GetOutput());
 
  const InternalPixelType lowerThreshold = std::stod("3600");
  const InternalPixelType upperThreshold = std::stod("4000");
  connectedThreshold->SetLower(lowerThreshold);
  connectedThreshold->SetUpper(upperThreshold);
  connectedThreshold->SetReplaceValue(255);
  InternalImageType::IndexType index1;
  index1[0] = std::stoi("295");
  index1[1] = std::stoi("98");
  index1[2] = std::stoi("159");
  connectedThreshold->AddSeed(index1);
 
  // InternalImageType::IndexType index2;
  // index2[0] = std::stoi("149");
  // index2[1] = std::stoi("172");
  // index2[2] = std::stoi("159");
  // connectedThreshold->AddSeed(index2);
 
  connectedThreshold->Update();
  //caster->SetInput(connectedThreshold->GetOutput());
  // writer->SetInput(connectedThreshold->GetOutput());
  // try
  // {
  //   writer->Update();
  // }
  // catch (const itk::ExceptionObject &excep)
  // {
  //   std::cerr << "Exception caught !" << std::endl;
  //   //std::cerr << excep << std::endl;
  // }
  cout << "保存完毕" << imgSize << endl;
  /
   //itk to vtk
    typedef itk::ImageToVTKImageFilter<InternalImageType> itkTovtkFilterType;
    itkTovtkFilterType::Pointer itkTovtkImageFilter = itkTovtkFilterType::New();
    itkTovtkImageFilter->SetInput(connectedThreshold->GetOutput());
    itkTovtkImageFilter->Update();
 
 
 
 
 vtkSmartPointer<vtkGPUVolumeRayCastMapper> volumeMapper =
        vtkSmartPointer<vtkGPUVolumeRayCastMapper>::New();
    volumeMapper->SetInputData(itkTovtkImageFilter->GetOutput());
    volumeMapper->SetSampleDistance(volumeMapper->GetSampleDistance() / 2);	//设置光线采样距离
 
    vtkSmartPointer<vtkVolumeProperty> volumeProperty =
        vtkSmartPointer<vtkVolumeProperty>::New();
    volumeProperty->SetInterpolationTypeToLinear();
    volumeProperty->ShadeOn();  //打开或者关闭阴影测试
    volumeProperty->SetAmbient(0.1);
    volumeProperty->SetDiffuse(0.9);  //漫反射
    volumeProperty->SetSpecular(0.2); //镜面反射
    volumeProperty->SetSpecularPower(10.0);
    //设置不透明度
    vtkSmartPointer<vtkPiecewiseFunction> opacityFun =
        vtkSmartPointer<vtkPiecewiseFunction>::New();
 
    opacityFun->AddPoint(-3024, 0);
    opacityFun->AddPoint(129.54, 0);
    opacityFun->AddPoint(145.24, 0.17);
    opacityFun->AddPoint(169.92, 0.63);
    opacityFun->AddPoint(395.58, 0.81);
    opacityFun->AddPoint(1578.73, 0.81);
    opacityFun->AddPoint(3071, 0.81);
 
 
    volumeProperty->SetScalarOpacity(opacityFun); //设置不透明度传输函数
 
    //设置颜色属性
    vtkSmartPointer<vtkColorTransferFunction> colorFun =
        vtkSmartPointer<vtkColorTransferFunction>::New();
    colorFun->AddRGBPoint(-3024, 0, 0, 0);
    colorFun->AddRGBPoint(129.54, 0.55, 0.25, 0.15);
    colorFun->AddRGBPoint(157.02, 1, 1, 1);
    colorFun->AddRGBPoint(169.92, 0.99, 0.87, 0.39);
    colorFun->AddRGBPoint(395.58, 1, 0.88, 0.66);
    colorFun->AddRGBPoint(1578.73, 1, 0.95, 0.96);
    colorFun->AddRGBPoint(3071, 0.83, 0.66, 1);
 
 
    volumeProperty->SetColor(colorFun);
    vtkSmartPointer<vtkVolume> volume =
        vtkSmartPointer<vtkVolume>::New();
    volume->SetMapper(volumeMapper);
    volume->SetProperty(volumeProperty);
 
    vtkSmartPointer<vtkCamera>camera = vtkSmartPointer<vtkCamera>::New();
    camera->SetViewUp(0, 0, -1);
    camera->SetPosition(0, 1, 0);
    camera->SetFocalPoint(0, 0, 0);
    camera->ComputeViewPlaneNormal();
    camera->Azimuth(10.0);
    camera->Elevation(10.0);
 
    // A renderer and render window
    vtkSmartPointer<vtkRenderer> renderer = vtkSmartPointer<vtkRenderer>::New();
    vtkSmartPointer<vtkRenderWindow> renderWindow = vtkSmartPointer<vtkRenderWindow>::New();
    renderWindow->AddRenderer(renderer);
    renderWindow->SetSize(1000, 1000);
    renderer->SetActiveCamera(camera);
    renderer->SetBackground(0, 0, 0);
    renderer->AddVolume(volume);
 
    // An interactor
    vtkSmartPointer<vtkRenderWindowInteractor> renderWindowInteractor =
        vtkSmartPointer<vtkRenderWindowInteractor>::New();
    renderWindowInteractor->SetRenderWindow(renderWindow);
    renderer->ResetCamera();
    renderWindow->Render();
    vtkSmartPointer<vtkInteractorStyleTrackballCamera>style = vtkSmartPointer<vtkInteractorStyleTrackballCamera>::New();
    renderWindowInteractor->SetInteractorStyle(style);
   // renderWindowInteractor->Initialize();
   // renderWindowInteractor->Start();
    renderWindowInteractor->Initialize();
    renderWindow->Render();
    renderWindowInteractor->Start();
 
  return EXIT_SUCCESS;
}
复制代码

重建出来的种植体