ITK 基于C++和ITK获取3D医学图像的最大连通区域

1、代码

#include "itkConnectedComponentImageFilter.h"
#include "itkImage.h"
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"
#include "itkImageRegionIterator.h"
#include "itkLabelShapeKeepNObjectsImageFilter.h"

using namespace std;

void print(const std::string& file, itk::Image<unsigned char, 2>* image)
{
    ofstream fout(file);
    auto size = image->GetBufferedRegion().GetSize();
    for (auto i = 0; i < size[0]; i++)
    {
        for (auto j = 0; j < size[1]; j++)
        {
            itk::Image<unsigned char, 2>::IndexType index{ {i,j} };
            auto c = image->GetPixel(index);
            fout << (int)c;
        }
        fout << std::endl;
    }
}


int main(int argc, char *argv[])
{
    // 1、获取二值化的分割结果mask
    const char* input_path = "D:/documents/vs2019/itk_tutorial/largestcompent/build/RelWithDebInfo/left_femur.nrrd";// 图像保存的路径
    using ImageType = itk::Image<unsigned char, 3>;
    using ReaderFilter = itk::ImageFileReader<ImageType>;
    auto reader = ReaderFilter::New();
    reader->SetFileName(input_path);
    reader->Update();
    ImageType::Pointer image = reader->GetOutput();

    // 2, 调用ITK的itkConnectedComponentImageFilter
    using ConnectedComponentFilter = itk::ConnectedComponentImageFilter<ImageType, ImageType>;
    auto connnectedComponentFilter = ConnectedComponentFilter::New();
    connnectedComponentFilter->SetInput(image);
    connnectedComponentFilter->Update();

    // 3,调用ITK的itkLabelShapeKeepNObjectsImageFilter
    using LabelShapeKeepNObjectsFilter = itk::LabelShapeKeepNObjectsImageFilter<ImageType>;
    auto labelShapeKeepNObjectsFilter = LabelShapeKeepNObjectsFilter::New();
    labelShapeKeepNObjectsFilter->SetInput(connnectedComponentFilter->GetOutput());
    labelShapeKeepNObjectsFilter->SetBackgroundValue(0);
    labelShapeKeepNObjectsFilter->SetNumberOfObjects(1);
    labelShapeKeepNObjectsFilter->SetAttribute(
        LabelShapeKeepNObjectsFilter::LabelObjectType::NUMBER_OF_PIXELS);
    labelShapeKeepNObjectsFilter->Update();
// 4、将Label便签转换为1
 typedef itk::RescaleIntensityImageFilter<ImageType, ImageType> RescaleFilterType;
  RescaleFilterType::Pointer rescaleFilter = RescaleFilterType::New();
  rescaleFilter->SetOutputMinimum(0);
  rescaleFilter->SetOutputMaximum(1);
  rescaleFilter->SetInput(labelShapeKeepNObjectsFilter->GetOutput());
  rescaleFilter->Update();
  //ImageType::Pointer largestccimage = rescaleFilter->GetOutput(); 

    // 5,保存最大连通区域图像

    const char* out_path = "D:/documents/vs2019/itk_tutorial/largestcompent/build/RelWithDebInfo/output.nii.gz"; // 保存输出结果的路径
    using WriterFilter = itk::ImageFileWriter<ImageType>;
    auto writer = WriterFilter::New();
    writer->SetFileName(out_path);
    writer->SetInput(rescaleFilter->GetOutput());
    writer->Update();
     
    return EXIT_SUCCESS;
}

2、效果

原始图像：

 去除小连通域之后

3、主要函数介绍

itk::ConnectedComponentImageFilter<ImageType, OutputImageType>;

该filter值检测连通区域，不编号，不排序。

itk::LabelShapeKeepNObjectsImageFilter;

该filter对连通区域进行排序，然后只保留最大的几个连通区域。

itk::LabelImageToShapeLabelMapFilter< OutputImageType, LabelMapType> ;

为保留的每个连通区域打上标签（标号），编号顺序为找到不同区域的顺序.

虽然LabelShapeKeepNObjectsImageFilter调用SetAttribute 之后，按照attribute对连通区域大小进行了排序，但该排序并不反映在label上。即大的label并不一定对应大的size。反映在 SetNumberOfObjects(n) 设定保留区域之后，只保留排序后前n个区域的label。

itk::RescaleIntensityImageFilter<OutputImageType, ImageType>；

上一步图像只有标签，图像上没有直接的区别，本步对保留下来的不同连通区域进行不同的染色，以示区分。