ITK学习笔记（五） ITK计算质心

通过LabelGeometryImageFilter可以求Label Image各个区域的质心和BoundingBox。然而，如果图像的spacing不是1，它并不会考虑进去，因此得到的结果并未我们所愿的。针对这个问题，有人实现了相关的模块(itkOBBLabelMap)，可以解决这个问题，将其添加到源码，然后重新配置和编译即可。（原文链接https://blog.csdn.net/Q1302182594/article/details/50447802）

QList <QPoint3D> calculateRegionCentroid(vtkImageData *data)
{
    QList < QPoint3D > centroids; // QPoint3D是仿照QPoint(不支持浮点数)而定义的一个类，支持浮点数
 
    if ( !data )
       return centroids;
 
    double spacing[3];
    data->GetSpacing( spacing );
 
    itk::VTKImageToImageFilter<UCharImageType>::Pointer   vtk2itkConnector // 将VTK图像转换为ITK图像
            = itk::VTKImageToImageFilter< UCharImageType>::New();
 
    BinaryImageToLabelMapFilterType::Pointer binaryImageToLabelMapFilter
            = BinaryImageToLabelMapFilterType::New();
 
    LabelMapToLabelImageFilterType::Pointer  labelMapToLabelImageFilter
            = LabelMapToLabelImageFilterType::New();
 
    LabelGeometryImageFilterType::Pointer labelGeometryImageFilter
            = LabelGeometryImageFilterType::New();
 
    vtk2itkConnector->SetInput( data );
 
    binaryImageToLabelMapFilter->SetInput( vtk2itkConnector->GetOutput() );
 
    labelMapToLabelImageFilter->SetInput( binaryImageToLabelMapFilter->GetOutput());
 
    labelGeometryImageFilter->SetInput( labelMapToLabelImageFilter->GetOutput() );
 
    // These generate optional outputs.
    labelGeometryImageFilter->CalculatePixelIndicesOn();
    labelGeometryImageFilter->CalculateOrientedBoundingBoxOn();
    labelGeometryImageFilter->CalculateOrientedLabelRegionsOn();
    labelGeometryImageFilter->Update();
 
    LabelGeometryImageFilterType::LabelsType allLabels = labelGeometryImageFilter->GetLabels();
    LabelGeometryImageFilterType::LabelsType::iterator allLabelsIt;
 
    for( allLabelsIt = allLabels.begin(); allLabelsIt != allLabels.end(); allLabelsIt++ ) {
        LabelGeometryImageFilterType::LabelPixelType labelValue = *allLabelsIt;
        if ( (int)labelValue > 0 ) {
            double x = labelGeometryImageFilter->GetCentroid( labelValue )[0];
            double y = labelGeometryImageFilter->GetCentroid( labelValue )[1];
            centroids.append( QPoint3D( x * spacing[0], y * spacing[1], 0 ) );
        }
    }
 
    binaryImageToLabelMapFilter->SetInput( NULL );
    labelMapToLabelImageFilter->SetInput( NULL );
    labelGeometryImageFilter->SetInput( NULL );
    vtk2itkConnector->SetInput( NULL );
 
    return centroids;
}

github资源：https://github.com/blowekamp/itkOBBLabelMap

int main(int, char *[])
{
vtkSmartPointer<vtkDICOMImageReader> reader = vtkSmartPointer<vtkDICOMImageReader>::New();
reader->SetFileName("C:\\Users\\zgb\\Desktop\\1-130.dcm");

/*************设置窗宽窗位************/

vtkSmartPointer<vtkImageMapToWindowLevelColors> imageViewer = vtkSmartPointer<vtkImageMapToWindowLevelColors>::New();
imageViewer->SetInputConnection(reader->GetOutputPort());
imageViewer->SetWindow(2894);
imageViewer->SetLevel(-81);
imageViewer->Update();

/*************** 阈值处理**************/

vtkSmartPointer<vtkImageThreshold> thresholdFilter = vtkSmartPointer<vtkImageThreshold>::New();
thresholdFilter->SetInputConnection(imageViewer->GetOutputPort());
thresholdFilter->ThresholdByUpper(200);
thresholdFilter->SetInValue(255);
thresholdFilter->SetOutValue(0);
thresholdFilter->Update();

vtkSmartPointer<vtkImageLuminance> luminanceFilter = vtkSmartPointer<vtkImageLuminance>::New();
luminanceFilter->SetInputConnection(thresholdFilter->GetOutputPort());
luminanceFilter->Update();

typedef itk::Image<unsigned char, 2> ImageTypetrans;
typedef itk::VTKImageToImageFilter<ImageTypetrans> VTKImageToImageType;

/************* vtk类型转换为itk类型 ************/

VTKImageToImageType::Pointer vtkImageToImageFilter = VTKImageToImageType::New();
vtkImageToImageFilter->SetInput(luminanceFilter->GetOutput());
vtkImageToImageFilter->Update();

typedef itk::Image<unsigned char, 2> ImageType;
ImageType::Pointer image = ImageType::New();


typedef itk::BinaryImageToLabelMapFilter<ImageType>
BinaryImageToLabelMapFilterType; //在一幅二值图像中标记连接组件，并产生一个标记组件的集合
BinaryImageToLabelMapFilterType::Pointer binaryImageToLabelMapFilter = BinaryImageToLabelMapFilterType::New();//实例化
binaryImageToLabelMapFilter->SetInput(vtkImageToImageFilter->GetOutput()); //输入图片
binaryImageToLabelMapFilter->Update(); //执行标记连接

typedef itk::LabelMapToLabelImageFilter<BinaryImageToLabelMapFilterType::OutputImageType, ImageType>
LabelMapToLabelImageFilterType; //将这些标记转化为标记的图像
LabelMapToLabelImageFilterType::Pointer labelMapToLabelImageFilter = LabelMapToLabelImageFilterType::New();//实例化
labelMapToLabelImageFilter->SetInput(binaryImageToLabelMapFilter->GetOutput());//输出作输入
labelMapToLabelImageFilter->Update(); //开始执行

typedef itk::LabelGeometryImageFilter< ImageType > LabelGeometryImageFilterType;
LabelGeometryImageFilterType::Pointer labelGeometryImageFilter = LabelGeometryImageFilterType::New();//实例化
labelGeometryImageFilter->SetInput( labelMapToLabelImageFilter->GetOutput() );//输入标记的图像
labelGeometryImageFilter->Update(); //开始执行

LabelGeometryImageFilterType::LabelsType allLabels =labelGeometryImageFilter->GetLabels(); //获取标记
LabelGeometryImageFilterType::LabelsType::iterator allLabelsIt;//定义迭代器

for( allLabelsIt=allLabels.begin();allLabelsIt!= allLabels.end(); allLabelsIt++ )

{
LabelGeometryImageFilterType::LabelPixelType labelValue = *allLabelsIt;//定义标记个数

std::cout << "\tCentroid: "
<< labelGeometryImageFilter->GetCentroid(labelValue) << std::endl;
std::cout << "Number of labels: "
<< labelGeometryImageFilter->GetNumberOfLabels() << std::endl;

}


return EXIT_SUCCESS;
}