Generate Structuring Elements With Accurate Area#
Synopsis#
Generate structuring elements with accurate area.
Results#
Output:
2D ball of radius 5 with radiusIsParametric mode off:
0 0 0 1 1 1 1 1 0 0 0
0 0 1 1 1 1 1 1 1 0 0
0 1 1 1 1 1 1 1 1 1 0
1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1
0 1 1 1 1 1 1 1 1 1 0
0 0 1 1 1 1 1 1 1 0 0
0 0 0 1 1 1 1 1 0 0 0
Expected foreground area: 78.5398
Computed foreground area: 97
Foreground area error: 23.5042%
2D ball of radius 5 with radiusIsParametric mode on:
0 0 0 0 0 1 0 0 0 0 0
0 0 1 1 1 1 1 1 1 0 0
0 1 1 1 1 1 1 1 1 1 0
0 1 1 1 1 1 1 1 1 1 0
0 1 1 1 1 1 1 1 1 1 0
1 1 1 1 1 1 1 1 1 1 1
0 1 1 1 1 1 1 1 1 1 0
0 1 1 1 1 1 1 1 1 1 0
0 1 1 1 1 1 1 1 1 1 0
0 0 1 1 1 1 1 1 1 0 0
0 0 0 0 0 1 0 0 0 0 0
Expected foreground area: 78.5398
Computed foreground area: 81
Foreground area error: 3.1324%
2D annulus of radius 5 and thickness 2 with radiusIsParametric mode off:
0 0 0 1 1 1 1 1 0 0 0
0 0 1 1 1 1 1 1 1 0 0
0 1 1 1 0 0 0 1 1 1 0
1 1 1 0 0 0 0 0 1 1 1
1 1 0 0 0 0 0 0 0 1 1
1 1 0 0 0 0 0 0 0 1 1
1 1 0 0 0 0 0 0 0 1 1
1 1 1 0 0 0 0 0 1 1 1
0 1 1 1 0 0 0 1 1 1 0
0 0 1 1 1 1 1 1 1 0 0
0 0 0 1 1 1 1 1 0 0 0
Expected foreground area: 50.2655
Computed foreground area: 60
Foreground area error: 19.3662%
2D annulus of radius 5 and thickness 2 with radiusIsParametric mode on:
0 0 0 0 0 1 0 0 0 0 0
0 0 1 1 1 1 1 1 1 0 0
0 1 1 1 1 0 1 1 1 1 0
0 1 1 0 0 0 0 0 1 1 0
0 1 1 0 0 0 0 0 1 1 0
1 1 0 0 0 0 0 0 0 1 1
0 1 1 0 0 0 0 0 1 1 0
0 1 1 0 0 0 0 0 1 1 0
0 1 1 1 1 0 1 1 1 1 0
0 0 1 1 1 1 1 1 1 0 0
0 0 0 0 0 1 0 0 0 0 0
Expected foreground area: 50.2655
Computed foreground area: 52
Foreground area error: 3.45071%
3D ball of radius 5 with radiusIsParametric mode off:
Expected foreground area: 523.599
Computed foreground area: 739
Foreground area error: 41.1386%
3D ball of radius 5 with radiusIsParametric mode on:
Expected foreground area: 523.599
Computed foreground area: 515
Foreground area error: 1.64224%
3D annulus of radius 5 and thickness 2 with radiusIsParametric mode off:
Expected foreground area: 410.501
Computed foreground area: 560
Foreground area error: 36.4185%
3D annulus of radius 5 and thickness 2 with radiusIsParametric mode on:
Expected foreground area: 410.501
Computed foreground area: 392
Foreground area error: 4.50703%
4D ball of radius 5 with radiusIsParametric mode off:
Expected foreground area: 3084.25
Computed foreground area: 4785
Foreground area error: 55.143%
4D ball of radius 5 with radiusIsParametric mode on:
Expected foreground area: 3084.25
Computed foreground area: 2929
Foreground area error: 5.03368%
4D annulus of radius 5 and thickness 2 with radiusIsParametric mode off:
Expected foreground area: 2684.53
Computed foreground area: 4024
Foreground area error: 49.8957%
4D annulus of radius 5 and thickness 2 with radiusIsParametric mode on:
Expected foreground area: 2684.53
Computed foreground area: 2504
Foreground area error: 6.72491%
Code#
C++#
#include "itkFlatStructuringElement.h"
// Helper function
template <class SEType>
bool
ComputeAreaError(SEType k, unsigned int thickness = 0);
int
main()
{
int scalarRadius = 5;
int scalarThickness = 2;
bool radiusIsParametric = true;
using SE2Type = itk::FlatStructuringElement<2>;
SE2Type::RadiusType r2;
r2.Fill(scalarRadius);
SE2Type k2;
std::cout << "2D ball of radius " << scalarRadius << " with radiusIsParametric mode off:" << std::endl;
k2 = SE2Type::Ball(r2);
ComputeAreaError(k2);
// Test the radiusIsParametric mode.
std::cout << "2D ball of radius " << scalarRadius << " with radiusIsParametric mode on:" << std::endl;
k2 = SE2Type::Ball(r2, radiusIsParametric);
ComputeAreaError(k2);
std::cout << "2D annulus of radius " << scalarRadius << " and thickness " << scalarThickness
<< " with radiusIsParametric mode off:" << std::endl;
k2 = SE2Type::Annulus(r2, scalarThickness, false);
ComputeAreaError(k2, scalarThickness);
// Test the radiusIsParametric mode.
std::cout << "2D annulus of radius " << scalarRadius << " and thickness " << scalarThickness
<< " with radiusIsParametric mode on:" << std::endl;
k2 = SE2Type::Annulus(r2, scalarThickness, false, radiusIsParametric);
ComputeAreaError(k2, scalarThickness);
using SE3Type = itk::FlatStructuringElement<3>;
SE3Type::RadiusType r3;
r3.Fill(scalarRadius);
SE3Type k3;
std::cout << "3D ball of radius " << scalarRadius << " with radiusIsParametric mode off:" << std::endl;
k3 = SE3Type::Ball(r3);
ComputeAreaError(k3);
// Test the radiusIsParametric mode.
std::cout << "3D ball of radius " << scalarRadius << " with radiusIsParametric mode on:" << std::endl;
k3 = SE3Type::Ball(r3, radiusIsParametric);
ComputeAreaError(k3);
std::cout << "3D annulus of radius " << scalarRadius << " and thickness " << scalarThickness
<< " with radiusIsParametric mode off:" << std::endl;
k3 = SE3Type::Annulus(r3, scalarThickness, false);
ComputeAreaError(k3, scalarThickness);
// Test the radiusIsParametric mode.
std::cout << "3D annulus of radius " << scalarRadius << " and thickness " << scalarThickness
<< " with radiusIsParametric mode on:" << std::endl;
k3 = SE3Type::Annulus(r3, scalarThickness, false, radiusIsParametric);
ComputeAreaError(k3, scalarThickness);
using SE4Type = itk::FlatStructuringElement<4>;
SE4Type::RadiusType r4;
r4.Fill(scalarRadius);
SE4Type k4;
std::cout << "4D ball of radius " << scalarRadius << " with radiusIsParametric mode off:" << std::endl;
k4 = SE4Type::Ball(r4);
ComputeAreaError(k4);
// Test the radiusIsParametric mode.
std::cout << "4D ball of radius " << scalarRadius << " with radiusIsParametric mode on:" << std::endl;
k4 = SE4Type::Ball(r4, radiusIsParametric);
ComputeAreaError(k4);
std::cout << "4D annulus of radius " << scalarRadius << " and thickness " << scalarThickness
<< " with radiusIsParametric mode off:" << std::endl;
k4 = SE4Type::Annulus(r4, scalarThickness, false);
ComputeAreaError(k4, scalarThickness);
// Test the radiusIsParametric mode.
std::cout << "4D annulus of radius " << scalarRadius << " and thickness " << scalarThickness
<< " with radiusIsParametric mode on:" << std::endl;
k4 = SE4Type::Annulus(r4, scalarThickness, false, radiusIsParametric);
ComputeAreaError(k4, scalarThickness);
return EXIT_SUCCESS;
}
template <class SEType>
bool
ComputeAreaError(SEType k, unsigned int thickness)
{
float expectedOuterForegroundArea = 1;
float expectedInnerForegroundArea;
if (thickness == 0)
{
// Circle/Ellipse has no inner area to subract.
expectedInnerForegroundArea = 0;
}
else
{
// Annulus does have inner area to subract.
expectedInnerForegroundArea = 1;
}
if (SEType::NeighborhoodDimension == 2)
{
expectedOuterForegroundArea *= itk::Math::pi;
expectedInnerForegroundArea *= itk::Math::pi;
}
else if (SEType::NeighborhoodDimension == 3)
{
expectedOuterForegroundArea *= 4.0 / 3.0 * itk::Math::pi;
expectedInnerForegroundArea *= 4.0 / 3.0 * itk::Math::pi;
}
else if (SEType::NeighborhoodDimension == 4)
{
expectedOuterForegroundArea *= 0.5 * itk::Math::pi * itk::Math::pi;
expectedInnerForegroundArea *= 0.5 * itk::Math::pi * itk::Math::pi;
}
else
{
return EXIT_FAILURE;
}
for (unsigned int i = 0; i < SEType::NeighborhoodDimension; ++i)
{
expectedOuterForegroundArea *= k.GetRadius()[i];
expectedInnerForegroundArea *= (k.GetRadius()[i] - thickness);
}
float expectedForegroundArea = expectedOuterForegroundArea - expectedInnerForegroundArea;
// Show the neighborhood if it is 2D.
typename SEType::Iterator SEIt;
if (SEType::NeighborhoodDimension == 2)
{
for (SEIt = k.Begin(); SEIt != k.End(); ++SEIt)
{
std::cout << *SEIt << "\t";
if ((SEIt - k.Begin() + 1) % k.GetSize()[0] == 0)
{
std::cout << std::endl;
}
}
}
// Compute the area/volume.
float computedForegroundArea = 0;
for (SEIt = k.Begin(); SEIt != k.End(); ++SEIt)
{
if (*SEIt)
{
computedForegroundArea++;
}
}
std::cout << "Expected foreground area: " << expectedForegroundArea << std::endl;
std::cout << "Computed foreground area: " << computedForegroundArea << std::endl;
std::cout << "Foreground area error: "
<< 100 * itk::Math::abs(expectedForegroundArea - computedForegroundArea) / expectedForegroundArea << "%"
<< "\n\n";
return EXIT_FAILURE;
}
Classes demonstrated#
-
template<unsigned int VDimension>
class FlatStructuringElement : public itk::Neighborhood<bool, VDimension> A class to support a variety of flat structuring elements, including versions created by decomposition of lines.
FlatStructuringElement provides several static methods, which can be used to create a structuring element with a particular shape, size, etc. Currently, those methods enable the creation of the following structuring elements: ball, box, cross, annulus, or polygon. Polygons are available as fast approximations of balls using line decompositions. Boxes also use line decompositions.
“Flat” refers to binary as opposed to grayscale structuring elements. Flat structuring elements can be used for both binary and grayscale images.
A Neighborhood has an N-dimensional radius. The radius is defined separately for each dimension as the number of pixels that the neighborhood extends outward from the center pixel. For example, a 2D Neighborhood object with a radius of 2x3 has sides of length 5x7. However, in the case of balls and annuli, this definition is slightly different from the parametric definition of those objects. For example, an ellipse of radius 2x3 has a diameter of 4x6, not 5x7. To have a diameter of 5x7, the radius would need to increase by 0.5 in each dimension. Thus, the “radius” of the neighborhood and the “radius” of the object should be distinguished.
To accomplish this, the “ball” and “annulus” structuring elements have an optional flag called “radiusIsParametric” (off by default). Setting this flag to true will use the parametric definition of the object and will generate structuring elements with more accurate areas, which can be especially important when morphological operations are intended to remove or retain objects of particular sizes. When the mode is turned off (default), the radius is the same, but the object diameter is set to (radius*2)+1, which is the size of the neighborhood region. Thus, the original ball and annulus structuring elements have a systematic bias in the radius of +0.5 voxels in each dimension relative to the parametric definition of the radius. Thus, we recommend turning this mode on for more accurate structuring elements, but this mode is turned off by default for backward compatibility.
As an example, a 3D ball of radius 5 should have an area of 523. With this mode turned on, the number of “on” pixels is 515 (error 1.6%), but with it turned off, the area is 739 (error 41%). For a 3D annulus of radius 5 and thickness 2, the area should be 410. With this mode turned on, the area is 392 (error 4.5%), but when turned off it is 560 (error 36%). This same trend holds for balls and annuli of any radius or dimension. For more detailed experiments with this mode, please refer to the results of the test itkFlatStructuringElementTest.cxx or the wiki example.
- ITK Sphinx Examples:
- ITK Sphinx Examples:
