ITK连通域分析

连通域

一个区域D中，任一闭曲线所围区域，都属于区域D，称为“单连通区域”。

连通域形状分析

1. 创建二值图像

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

ImageType::Pointer CreateMaskImage()
{
	ImageType::Pointer image = ImageType::New();
	ImageType::IndexType start;
	start.Fill(0);

	ImageType::SizeType size;
	size.Fill(128);

	ImageType::RegionType region;
	region.SetSize(size);
	region.SetIndex(start);

	image->SetRegions(region);
	image->Allocate();

	using IteratorType = itk::ImageRegionIteratorWithIndex< ImageType >;
	IteratorType it(image, image->GetRequestedRegion());

	for (it.GoToBegin(); !it.IsAtEnd(); ++it)
	{
		ImageType::IndexType idx = it.GetIndex();

		// 创建空心圆 C=（50,50）,R=15,r=5
		if ((idx[0] - 50)*(idx[0] - 50) + (idx[1] - 50)*(idx[1] - 50) < 225 &&
			(idx[0] - 50)*(idx[0] - 50) + (idx[1] - 50)*(idx[1] - 50) > 25)
			it.Set(1);

		// 创建圆 C=(70,80),r=7
		else if ((idx[0] - 70)*(idx[0] - 70) + (idx[1] - 80)*(idx[1] - 80) < 49)
			it.Set(1);

		// 创建矩形 C=(100,100),W(120,110)
		else if (idx[0] >= 100 && idx[0] <= 120 &&
			idx[1] >= 100 && idx[1] <= 110)
			it.Set(1);

		else
			it.Set(0);
	}

	return image;
}

生成图像如下

2. 打印所有连通区域信息

typedef unsigned long LabelType;
typedef itk::ShapeLabelObject< LabelType, 2 >									LabelObjectType;
typedef itk::LabelMap< LabelObjectType >										LabelMapType;
typedef itk::BinaryImageToShapeLabelMapFilter< ImageType, LabelMapType >		I2LType;
I2LType::Pointer i2l = I2LType::New();
i2l->SetInput(mask);
i2l->SetFullyConnected(false);
i2l->SetInputForegroundValue(1);
i2l->SetOutputBackgroundValue(0);
i2l->Update();

LabelMapType::Pointer labelMap = i2l->GetOutput();
for (unsigned int label = 1; label <= labelMap->GetNumberOfLabelObjects(); label++)
{
	const LabelObjectType * labelObject = labelMap->GetLabelObject(label);
	labelObject->Print(std::cout);
	std::cout << "----------------------------------------" << std::endl;
}

打印信息如下：

ShapeLabelObject (000002469EAA6CE0)
  RTTI typeinfo:   class itk::ShapeLabelObject<unsigned long,2>
  Reference Count: 1
  LineContainer: 000002469EAA6CF0
  Label: 1
  NumberOfPixels: 616
  PhysicalSize: 616
  Perimeter: 127.254
  NumberOfPixelsOnBorder: 0
  PerimeterOnBorder: 0
  PerimeterOnBorderRatio: 0
  Elongation: 1
  Flatness: 1
  Roundness: 0.691393
  Centroid: [50, 50]
  BoundingBox:   ImageRegion (000002469EAA6D20)
    Dimension: 2
    Index: [36, 36]
    Size: [29, 29]
  EquivalentSphericalRadius: 14.0028
  EquivalentSphericalPerimeter: 87.9823
  EquivalentEllipsoidDiameter: [28.0056, 28.0056]
  PrincipalMoments: [61.9156, 61.9156]
  PrincipalAxes:
1 0
0 1
  FeretDiameter: 0
  m_OrientedBoundingBoxSize: [0, 0]
  m_OrientedBoundingBoxOrigin: [0, 0]
----------------------------------------
ShapeLabelObject (000002469EAA7030)
  RTTI typeinfo:   class itk::ShapeLabelObject<unsigned long,2>
  Reference Count: 1
  LineContainer: 000002469EAA7040
  Label: 2
  NumberOfPixels: 145
  PhysicalSize: 145
  Perimeter: 41.524
  NumberOfPixelsOnBorder: 0
  PerimeterOnBorder: 0
  PerimeterOnBorderRatio: 0
  Elongation: 1
  Flatness: 1
  Roundness: 1.02799
  Centroid: [70, 80]
  BoundingBox:   ImageRegion (000002469EAA7070)
    Dimension: 2
    Index: [64, 74]
    Size: [13, 13]
  EquivalentSphericalRadius: 6.79374
  EquivalentSphericalPerimeter: 42.6863
  EquivalentEllipsoidDiameter: [13.5875, 13.5875]
  PrincipalMoments: [11.5172, 11.5172]
  PrincipalAxes:
1 0
0 1
  FeretDiameter: 0
  m_OrientedBoundingBoxSize: [0, 0]
  m_OrientedBoundingBoxOrigin: [0, 0]
----------------------------------------
ShapeLabelObject (000002469EAA7170)
  RTTI typeinfo:   class itk::ShapeLabelObject<unsigned long,2>
  Reference Count: 1
  LineContainer: 000002469EAA7180
  Label: 3
  NumberOfPixels: 231
  PhysicalSize: 231
  Perimeter: 59.5651
  NumberOfPixelsOnBorder: 0
  PerimeterOnBorder: 0
  PerimeterOnBorderRatio: 0
  Elongation: 1.91485
  Flatness: 1.91485
  Roundness: 0.904522
  Centroid: [110, 105]
  BoundingBox:   ImageRegion (000002469EAA71B0)
    Dimension: 2
    Index: [100, 100]
    Size: [21, 11]
  EquivalentSphericalRadius: 8.57494
  EquivalentSphericalPerimeter: 53.8779
  EquivalentEllipsoidDiameter: [12.3935, 23.7317]
  PrincipalMoments: [10, 36.6667]
  PrincipalAxes:
0 1
-1 -0
  FeretDiameter: 0
  m_OrientedBoundingBoxSize: [0, 0]
  m_OrientedBoundingBoxOrigin: [0, 0]
----------------------------------------

3. 统计灰度信息

typedef itk::BinaryImageToStatisticsLabelMapFilter< ImageType, ImageType > ConverterType;
ConverterType::Pointer converter = ConverterType::New();
converter->SetInput(reader->GetOutput());
converter->SetFeatureImage(reader2->GetOutput());
converter->SetInputForegroundValue(1);
converter->SetFullyConnected(false);
converter->Update();
converter->GetOutput()->PrintLabelObjects();

4. 过滤连通区域

typedef unsigned long LabelType;
typedef itk::ShapeLabelObject< LabelType, 2 >									LabelObjectType;
typedef itk::LabelMap< LabelObjectType >										LabelMapType;
typedef itk::BinaryImageToShapeLabelMapFilter< ImageType, LabelMapType >		I2LType;
I2LType::Pointer i2l = I2LType::New();
i2l->SetInput(mask);
i2l->SetFullyConnected(false);
i2l->SetInputForegroundValue(1);
i2l->SetOutputBackgroundValue(0);
i2l->Update();

LabelMapType::Pointer labelMap = i2l->GetOutput();
for (unsigned int label = 1; label <= labelMap->GetNumberOfLabelObjects(); label++)
{
	const LabelObjectType * labelObject = labelMap->GetLabelObject(label);
	labelObject->Print(std::cout);
	std::cout << "----------------------------------------" << std::endl;
}

std::vector< LabelObjectType::Pointer > filtered_label_objects;
for (LabelMapType::Iterator it(labelMap); !it.IsAtEnd(); ++it)
{
	// 过滤条件
	if (it.GetLabelObject()->GetRoundness() < 0.95)
		filtered_label_objects.push_back(it.GetLabelObject());
}

for (int i = 0; i < filtered_label_objects.size(); i++)
{
	labelMap->RemoveLabelObject(filtered_label_objects[i]);
}

typedef itk::LabelMapToBinaryImageFilter< LabelMapType, ImageType> L2IType;
L2IType::Pointer l2i = L2IType::New();
l2i->SetInput(labelMap);
l2i->Update();

using WriterType = itk::ImageFileWriter<ImageType>;
WriterType::Pointer writer = WriterType::New();
writer->SetFileName("mask.mha");
writer->SetInput(l2i->GetOutput());
writer->SetImageIO(itk::MetaImageIO::New());
writer->Update();

效果如下：

总结

通过连通域形状分析可以过滤到阈值分割中的错误部分，另外ITK还提供了itk::BinaryImageToStatisticsLabelMapFilter类，可以提供原图区域的统计信息。

itk::ShapeLabelObject提供的形状信息包括：

1. NumberOfPixels 像素数
2. PhysicalSize 物理尺寸
3. Perimeter 周长
4. NumberOfPixelsOnBorder 边界像素数
5. PerimeterOnBorder 边界周长
6. PerimeterOnBorderRatio 边界周长比
7. Elongation 伸长率
8. Flatness 平整度
9. Roundness 圆度
10. Centroid 质心
11. BoundingBox 包围盒
12. EquivalentSphericalRadius 等效球面半径
13. EquivalentSphericalPerimeter 等效球面周长
14. EquivalentEllipsoidDiameter 等效椭球直径
15. PrincipalMoments 主力矩
16. PrincipalAxes 主轴
17. FeretDiameter

itk::StatisticsLabelObject提供的统计信息包括：

1. Minimum
2. Maximum
3. Mean
4. Sum
5. StandardDeviation 标准差
6. Variance 方差
7. Median
8. Skewness 偏态
9. Kurtosis 峰度
10. WeightedElongation 加权伸长率
11. WeightedFlatness 加权平面度
12. MaximumIndex
13. MinimumIndex
14. CenterOfGravity 重心
15. WeightedPrincipalMoments 加权主力矩
16. WeightedPrincipalAxes 加权主轴