ITK 合并 dicom image 图层合并,合并两个图层内容到新的图层

 要求：要求图像空间大小必须是一样的，要不没办法合并，不紧紧是 Spacing，Dimensions

原图层（这里我截取的大小不同，实际必须相同）

合并后

 还可以实现，异或，差；

/*=========================================================================
 *
 *  Copyright NumFOCUS
 *
 *  Licensed under the Apache License, Version 2.0 (the "License");
 *  you may not use this file except in compliance with the License.
 *  You may obtain a copy of the License at
 *
 *         https://www.apache.org/licenses/LICENSE-2.0.txt
 *
 *  Unless required by applicable law or agreed to in writing, software
 *  distributed under the License is distributed on an "AS IS" BASIS,
 *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 *  See the License for the specific language governing permissions and
 *  limitations under the License.
 *
 *=========================================================================*/
#include "itkImage.h"
#include "itkRescaleIntensityImageFilter.h"
#include "itkAddImageFilter.h"

#include <itkImageToVTKImageFilter.h>

#include "vtkVersion.h"
#include "vtkImageViewer.h"
#include "vtkImageMapper3D.h"
#include "vtkRenderWindowInteractor.h"
#include "vtkSmartPointer.h"
#include "vtkImageActor.h"
#include "vtkInteractorStyleImage.h"
#include "vtkRenderer.h"
#include "itkSymmetricSecondRankTensor.h"
#include "itkGDCMImageIO.h"
#include "itkImage.h"
#include "itkImageSeriesReader.h"
#include "itkGDCMSeriesFileNames.h"
#include "itkImageSeriesWriter.h"
#include "itkImageFileWriter.h"
#include "itkGDCMImageIOFactory.h"

#include "itkVTKImageIOFactory.h"

using namespace std;

string directoryName = "D:/datasource/baimingxue/femur-cartilage/ScalarVolume_25"; 
string directoryName1 = "D:/datasource/baimingxue/femur/ScalarVolume_17";
using PixelType = signed short;
constexpr unsigned int Dimension = 3;

//using ImageType = itk::Image<unsigned char, Dimension>;
 

int
main()
{
  //初始化待读取序列的格式类型
  using ImageType = itk::Image<PixelType, Dimension>;
  using ReaderType = itk::ImageSeriesReader<ImageType>;
  using ImageIOType = itk::GDCMImageIO;
  using NamesGeneratorType = itk::GDCMSeriesFileNames;

  //设置IO，并获取文件名
  ImageIOType::Pointer        gdcmIO = ImageIOType::New();
  NamesGeneratorType::Pointer namesGenerator = NamesGeneratorType::New();
  namesGenerator->SetInputDirectory(directoryName);
  const ReaderType::FileNamesContainer& filenames =
    namesGenerator->GetInputFileNames();

  //输出文件名
  std::size_t numberOfFileNames = filenames.size();
  std::cout << numberOfFileNames << std::endl;
  for (unsigned int fni = 0; fni < numberOfFileNames; ++fni)
  {
    std::cout << "filename # " << fni << " = ";
    std::cout << filenames[fni] << std::endl;
  }

  //读取dcm序列
  ReaderType::Pointer reader = ReaderType::New();
  reader->SetImageIO(gdcmIO);
  reader->SetFileNames(filenames);
  try
  {
    reader->Update();
  }
  catch (const itk::ExceptionObject& e)
  {
    std::cerr << "exception in file reader " << std::endl;
     
    return EXIT_FAILURE;
  }

  //读取dcm序列
  ReaderType::Pointer reader1 = ReaderType::New();

  {

    //设置IO，并获取文件名
    ImageIOType::Pointer        gdcmIO = ImageIOType::New();
    NamesGeneratorType::Pointer namesGenerator = NamesGeneratorType::New();
    namesGenerator->SetInputDirectory(directoryName);
    const ReaderType::FileNamesContainer& filenames =
      namesGenerator->GetInputFileNames();

    //输出文件名
    std::size_t numberOfFileNames = filenames.size();
    std::cout << numberOfFileNames << std::endl;
    for (unsigned int fni = 0; fni < numberOfFileNames; ++fni)
    {
      std::cout << "filename # " << fni << " = ";
      std::cout << filenames[fni] << std::endl;
    }

    //读取dcm序列
 
    reader1->SetImageIO(gdcmIO);
    reader1->SetFileNames(filenames);
    try
    {
      reader1->Update();
    }
    catch (const itk::ExceptionObject& e)
    {
      std::cerr << "exception in file reader " << std::endl;

      return EXIT_FAILURE;
    }
  }
 

  using AddImageFilterType = itk::AddImageFilter<ImageType, ImageType>;

  auto addFilter = AddImageFilterType::New();
  addFilter->SetInput1(reader->GetOutput());
  addFilter->SetInput2(reader1->GetOutput());
  addFilter->Update();

  // Visualize first image
  using ConnectorType = itk::ImageToVTKImageFilter<ImageType>;
  auto connector1 = ConnectorType::New();
  connector1->SetInput(reader->GetOutput());

  vtkSmartPointer<vtkImageActor> actor1 = vtkSmartPointer<vtkImageActor>::New();
#if VTK_MAJOR_VERSION <= 5
  actor1->SetInput(connector1->GetOutput());
#else
  connector1->Update();
  actor1->GetMapper()->SetInputData(connector1->GetOutput());
#endif
  // Visualize first image
  using ConnectorType = itk::ImageToVTKImageFilter<ImageType>;
  auto connector2 = ConnectorType::New();
  connector2->SetInput(reader1->GetOutput());

  vtkSmartPointer<vtkImageActor> actor2 = vtkSmartPointer<vtkImageActor>::New();
#if VTK_MAJOR_VERSION <= 5
  actor2->SetInput(connector2->GetOutput());
#else
  connector2->Update();
  actor2->GetMapper()->SetInputData(connector2->GetOutput());
#endif

  // write file 
  itk::VTKImageIOFactory::RegisterOneFactory();
  typedef itk::ImageFileWriter< ImageType > WriterType;
  WriterType::Pointer writer = WriterType::New();
  writer->SetFileName("result_3.mhd");
  writer->SetInput(addFilter->GetOutput());
  std::cout << "Writing the image as " << std::endl << std::endl;
 
  try
  {
    writer->Update();
  }
  catch (itk::ExceptionObject &ex)
  {
    std::cout << ex << std::endl;
    return EXIT_FAILURE;
  }
  // Visualize joined image
  auto addConnector = ConnectorType::New();
  addConnector->SetInput(addFilter->GetOutput());

  vtkSmartPointer<vtkImageActor> addActor = vtkSmartPointer<vtkImageActor>::New();
#if VTK_MAJOR_VERSION <= 5
  addActor->SetInput(addConnector->GetOutput());
#else
  addConnector->Update();
  addActor->GetMapper()->SetInputData(addConnector->GetOutput());
#endif
  // There will be one render window
  vtkSmartPointer<vtkRenderWindow> renderWindow = vtkSmartPointer<vtkRenderWindow>::New();
  renderWindow->SetSize(900, 300);

  vtkSmartPointer<vtkRenderWindowInteractor> interactor = vtkSmartPointer<vtkRenderWindowInteractor>::New();
  interactor->SetRenderWindow(renderWindow);

  // Define viewport ranges
  // (xmin, ymin, xmax, ymax)
  double leftViewport[4] = { 0.0, 0.0, 0.33, 1.0 };
  double centerViewport[4] = { 0.33, 0.0, 0.66, 1.0 };
  double rightViewport[4] = { 0.66, 0.0, 1.0, 1.0 };

  // Setup both renderers
  vtkSmartPointer<vtkRenderer> leftRenderer = vtkSmartPointer<vtkRenderer>::New();
  renderWindow->AddRenderer(leftRenderer);
  leftRenderer->SetViewport(leftViewport);
  leftRenderer->SetBackground(.6, .5, .4);

  vtkSmartPointer<vtkRenderer> centerRenderer = vtkSmartPointer<vtkRenderer>::New();
  renderWindow->AddRenderer(centerRenderer);
  centerRenderer->SetViewport(centerViewport);
  centerRenderer->SetBackground(.4, .5, .6);

  vtkSmartPointer<vtkRenderer> rightRenderer = vtkSmartPointer<vtkRenderer>::New();
  renderWindow->AddRenderer(rightRenderer);
  rightRenderer->SetViewport(rightViewport);
  rightRenderer->SetBackground(.4, .5, .6);

  // Add the sphere to the left and the cube to the right
  leftRenderer->AddActor(actor1);
  centerRenderer->AddActor(actor2);
  rightRenderer->AddActor(addActor);

  leftRenderer->ResetCamera();
  centerRenderer->ResetCamera();
  rightRenderer->ResetCamera();

  renderWindow->Render();

  vtkSmartPointer<vtkInteractorStyleImage> style = vtkSmartPointer<vtkInteractorStyleImage>::New();
  interactor->SetInteractorStyle(style);

  interactor->Start();

  return EXIT_SUCCESS;
}