How can I apply a 2D mask in a 3D image? - SliceBySliceImageFilter

Hi!

I am trying to mask a 3D image with a 2D mask with the itkSliceBySliceImageFilter. Here is my code:

// We create the filter mask
using MaskFilterType = itk::MaskImageFilter<ITKImgLabelType2D, ITKImgLabelType2D, ITKImgLabelType2D>;
MaskFilterType::Pointer maskFilter3D = MaskFilterType::New();
maskFilter3D->SetInput(image);
maskFilter3D->SetMaskImage(mask);
maskFilter3D->Update();

// Apply a 2D mask to every slice of a 3D image
typedef itk::SliceBySliceImageFilter < ITKImgInType, ITKImgLabelType, MaskFilterType, MaskFilterType, ITKImgLabelType2D, ITKImgLabelType2D > SliceBySliceImageFilter;
SliceBySliceImageFilter::Pointer sliceBySliceFilter = SliceBySliceImageFilter::New();
sliceBySliceFilter->SetFilter(maskFilter3D);
sliceBySliceFilter->SetInput(&image3D);

try
{
sliceBySliceFilter->Update();
}
catch (const itk::ExceptionObject& err)

Where 'ITKImgLabelType2D = itk::Image<unsigned char, 2>ā€™ , ā€˜ITKImgInType = itk::Image<short, 3>ā€™ and ā€˜ITKImgLabelType = itk::Image<unsigned char, 3>ā€™

When I execute sliceBySliceFilter->Update(); I get the following error:

ā€˜itk::ERROR: SliceBySliceImageFilter(000001630B273560): At least 2 of the first 2 indexed inputs are required but only 1 are specified. The required inputs are expected to be the first inputsā€™

And I have no idea of what could it be. How can I solve it? Could I obtain what I want in a different way?

My ITK version is: 4.8.2

Thank you very much,
Cristina

Maybe try sliceBySliceFilter->SetInput(1, mask);? Looking at the code, SliceBySliceImageFilter should be able to work with multiple inputs.

But it expects all the inputs to have the same size (requested region). Depending on what you are trying to accomplish, it might be easier to construct a 3D mask using TileImageFilter, and use 3D variant of MaskImageFilter.

Thank you for the answer.

  • Firstly

I have tried the first advice:

auto mask = ITKImgLabelType2D::New();
ITKImgLabelType2D::RegionType region = a_pmask->GetLargestPossibleRegion();

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

ITKImgLabelType2D::SizeType regionSize = region.GetSize();
itk::ImageRegionIterator imageIterator(mask, region);

while (!imageIterator.IsAtEnd())
{
imageIterator.Set(itk::NumericTraits::max());
++imageIterator;
}

// We create the filter mask
using MaskFilterType = itk::MaskImageFilter<ITKImgLabelType2D, ITKImgLabelType2D, ITKImgLabelType2D>;
MaskFilterType::Pointer maskFilter3D = MaskFilterType::New();
maskFilter3D->SetInput(a_pmask);
maskFilter3D->SetMaskImage(mask);

// Apply a 2D mask to every slice of a 3D image
typedef itk::SliceBySliceImageFilter< ITKImgInType, ITKImgLabelType, MaskFilterType, MaskFilterType, ITKImgLabelType2D, ITKImgLabelType2D > SliceBySliceImageFilter;
SliceBySliceImageFilter::Pointer sliceBySliceFilter = SliceBySliceImageFilter::New();
sliceBySliceFilter->SetFilter(maskFilter3D);
sliceBySliceFilter->SetInput(1, &a_pimage3D);

And I get another error, which is:

itk::ERROR: SliceBySliceImageFilter(0000019653A533C0): Input Primary is required but not set.


  • Secondly

I have tried to construct a 3D mask using TileImageFilter , and use 3D variant of MaskImageFilter.

Here is my code:

constexpr unsigned int InputDimension = 2;
constexpr unsigned int OutputDimension = 3;

using PixelType = unsigned char;
using InputImageType = itk::Image<PixelType, InputDimension>;
using OutputImageType = itk::Image<PixelType, OutputDimension>;

using FilterType = itk::TileImageFilter<InputImageType, OutputImageType>;
auto filter = FilterType::New();

itk::FixedArray<unsigned int, OutputDimension> layout;
layout[0] = 2;
layout[1] = 2;
layout[2] = 0;

filter->SetLayout(layout);

ITKImgInType::SizeType insize = a_pimage3D.GetRequestedRegion().GetSize();
int nImages = insize[2];

for (int ii = 0; ii <= nImages - 1; ++ii)
{
InputImageType::Pointer input = a_pmask;
filter->SetInput(ii, input);
}

constexpr PixelType defaultValue = 128;
filter->SetDefaultPixelValue(defaultValue);

filter->GetOutput()->SetSpacing(a_pimage3D.GetSpacing());
filter->GetOutput()->SetOrigin(a_pimage3D.GetOrigin());
filter->GetOutput()->SetDirection(a_pimage3D.GetDirection());

if (f_verbose)
{
a_file << "Spacing 3D mask after: " << filter->GetOutput()->GetSpacing() << std::endl;
a_file << "Spacing 3D image: " << a_pimage3D.GetSpacing() << std::endl << std::endl;
a_file << "Origin 3D mask after: " << filter->GetOutput()->GetOrigin() << std::endl;
a_file << "Origin 3D image: " << a_pimage3D.GetOrigin() << std::endl << std::endl;
}

// Apply a 3D mask to a 3D image
using MaskFilterType = itk::MaskImageFilter<ITKImgInType, ITKImgLabelType, ITKImgLabelType>;
MaskFilterType::Pointer maskFilter3D = MaskFilterType::New();
maskFilter3D->SetInput(&a_pimage3D);
maskFilter3D->SetMaskImage(filter->GetOutput());

try
{
maskFilter3D->Update();
}
catch (const itk::ExceptionObject& err)

I donā€™t understand whatā€™s wrong, because if I print the origin and the spacing of the a_image3D and the filter, as you can see in my code, I get the following information:

Spacing 3D mask after: [0.847656, 0.847656, 2.98663]
Spacing 3D image: [0.847656, 0.847656, 2.98663]

Origin 3D mask after: [-208.208, -251.595, 1148.21]
Origin 3D image: [-208.208, -251.595, 1148.21]

So, the spacing and the origin is the same. But then, when I execute maskFilter3D->Update(); I get the following error:

ExceptionObject caught !Alma::SegmentationLib::ApplyMask3D
C:\Dev\ITK\Modules\Core\Common\include\itkImageToImageFilter.hxx:248:
itk::ERROR: MaskImageFilter(00000231EBE04390): Inputs do not occupy the same physical space!
InputImage Origin: [-2.0820781e+02, -2.5159518e+02, 1.1482096e+03], InputImage_1 Origin: [-2.0820781e+02, -2.5159518e+02, 0.0000000e+00]
Tolerance: 8.4765559e-07
InputImage Spacing: [8.4765559e-01, 8.4765559e-01, 2.9866300e+00], InputImage_1 Spacing: [8.4765559e-01, 8.4765559e-01, 1.0000000e+00]
Tolerance: 8.4765559e-07

What should I do?

Thank you very much,
Cristina

From the docs: ā€œThe output image have a larger dimension than the input images. This filter can be used to create a volume from a series of inputs by specifying a layout of 1,1,0.ā€ You are using a tile layout 2,2,0. Your mask is probably bigger than the image.

a_file << "Spacing 3D mask after: " << filter->GetOutput()->GetSpacing() << std::endl;

You should print to the same or higher precision as the comparison tolerance (7-8 decimal places) in order to spot a possible problem here.

Thanks!

Iā€™ve changed the tile layout and I also get exactly the same error.

My code now is:

using PixelType = unsigned char;
using InputImageType = ITKImgLabelType2D;
using OutputImageType = ITKImgLabelType;

using FilterType = itk::TileImageFilter<InputImageType, OutputImageType>;
auto filter = FilterType::New();

itk::FixedArray<unsigned int, 3> layout;
layout[0] = 1;
layout[1] = 1;
layout[2] = 0;

filter->SetLayout(layout);

ITKImgInType::SizeType insize = a_pimage3D.GetRequestedRegion().GetSize();
int nImages = insize[2];

unsigned int inputImageNumber = 0;

for (int ii = 1; ii < nImages-1; ++ii)
{
InputImageType::Pointer input = a_pmask;
input->DisconnectPipeline();
filter->SetInput(inputImageNumber++, input);
}

constexpr PixelType defaultValue = 128;
filter->SetDefaultPixelValue(defaultValue);

filter->GetOutput()->SetSpacing(a_pimage3D.GetSpacing());
filter->GetOutput()->SetOrigin(a_pimage3D.GetOrigin());
filter->GetOutput()->SetDirection(a_pimage3D.GetDirection());

if (f_verbose)
{
a_file << std::fixed << std::setprecision(9) << std::endl;
a_file << "Spacing 3D mask after: " << filter->GetOutput()->GetSpacing() << std::endl;
a_file << "Spacing 3D image: " << a_pimage3D.GetSpacing() << std::endl << std::endl;
a_file << "Origin 3D mask after: " << filter->GetOutput()->GetOrigin() << std::endl;
a_file << "Origin 3D image: " << a_pimage3D.GetOrigin() << std::endl << std::endl;
}

// Apply a 3D mask to a 3D image
using MaskFilterType = itk::MaskImageFilter<ITKImgInType, ITKImgLabelType, ITKImgLabelType>;
MaskFilterType::Pointer maskFilter3D = MaskFilterType::New();
maskFilter3D->SetInput1(&a_pimage3D);
maskFilter3D->SetInput2(filter->GetOutput());

try
{
maskFilter3D->Update();
}
catch (const itk::ExceptionObject& err)

What exactly do the layout values ā€‹ā€‹1,1,0 mean?

Could it be possible that the error is related to the different input pixels type of the MaskImageFilter function?

Thank you,
Cristina

Have you tried filter->Update(); before if (f_verbose)? Are size, bufferedRegion and direction same, too? You could inspect that using a debugger, instead of printing to console.

Yes, I have tried that. And I also have checked that the inputs have the same size, bufferedRegion and direction too.

As a possible solution, I have changed the way of resampling the image.
I have replaced this:

filter->GetOutput()->SetSpacing(a_pimage3D.GetSpacing());
filter->GetOutput()->SetOrigin(a_pimage3D.GetOrigin());
filter->GetOutput()->SetDirection(a_pimage3D.GetDirection());

For this:

using ResampleFilterType = itk::ResampleImageFilter<ITKImgLabelType, ITKImgLabelType>;
auto resampleFilter = ResampleFilterType::New();

resampleFilter->SetInput(filter->GetOutput());
resampleFilter->SetOutputSpacing(a_pimage3D.GetSpacing());
resampleFilter->SetOutputOrigin(a_pimage3D.GetOrigin());
resampleFilter->SetOutputDirection(a_pimage3D.GetDirection());

This way, I donā€™t get the error: Inputs do not occupy the same physical space!, but I do get the following error: Requested region is (at least partially) outside the largest possible region.

Do you know the reason for this?

Thank you for your dedication.
Cristina

You need to invoke filter->Update(), or at least filter->UpdateOutputInformation() before filter->GetOutput() will have meaningful data.

Hi!

I have tried that and it didnā€™t work for me. But I have also tried another thing, that consists in ā€œresamplingā€ the mask with the ResampleImageFilter instead of manually. This way, I have fixed the error. But now another thing is happening.

I create the 3D mask correctly, but then, when I ā€œresampleā€ it, the 3D mask disappears.

Iā€™m going to explain myself with my code.

First of all I create a 3D mask from the 2D mask:

using PixelType = unsigned char;
using InputImageType = ITKImgLabelType2D;
using OutputImageType = ITKImgLabelType;

using FilterType = itk::TileImageFilter<InputImageType, OutputImageType>;
auto filter = FilterType::New();

itk::FixedArray<unsigned int, 3> layout;
layout[0] = 1;
layout[1] = 1;
layout[2] = 0;

filter->SetLayout(layout);

ITKImgInType::SizeType insize = a_pimage3D.GetRequestedRegion().GetSize();
int nImages = insize[2];

unsigned int inputImageNumber = 0;
InputImageType::Pointer input = a_pmask;

WriteImage(a_pmask, ā€œa_pmask.jpgā€, a_file);
for (int ii = 0; ii < nImages; ++ii)
{
input->DisconnectPipeline();
filter->SetInput(inputImageNumber++, a_pmask);
}

try
{
filter->UpdateOutputInformation();
filter->UpdateLargestPossibleRegion();
}
catch (const itk::ExceptionObject& err)
{
if (f_verbose)
{
a_file << ā€œExceptionObject caught !ā€ << func << std::endl;
a_file << err.what() << std::endl;
}
}

Thatā€™s done correctly, because if I print the image I can see the 2D mask in every slice of the 3D mask.
Then, I ā€œresampleā€ the mask:

ITKImgInType::SizeType outputSize = a_pimage3D.GetLargestPossibleRegion().GetSize();
ITKImgLabelType::SpacingType outputSpacing = a_pimage3D.GetSpacing();

using ScalarType = double;

using NearestNeighborInterpolatorType = itk::NearestNeighborInterpolateImageFunction<ITKImgLabelType, ScalarType>;
auto nearestNeighborInterpolator = NearestNeighborInterpolatorType::New();

using ResampleImageFilterType = itk::ResampleImageFilter<ITKImgLabelType, ITKImgLabelType>;
ResampleImageFilterType::Pointer resample = ResampleImageFilterType::New();

resample->SetInput(filter->GetOutput());
resample->SetReferenceImage(&a_pimage3D);
resample->SetInterpolator(nearestNeighborInterpolator);
resample->UseReferenceImageOn(); // spacing, origin and direction of the reference image will be used
resample->SetSize(outputSize);

try
{
resample->UpdateOutputInformation();
resample->UpdateLargestPossibleRegion();
}
catch (const itk::ExceptionObject& err)
{
if (f_verbose)
{
a_file << ā€œExceptionObject caught !ā€ << func << std::endl;
a_file << err.what() << std::endl;
}
}

Now, if I print the 3D mask, is wholly or partially black. So, when I do apply this mask to my original 3D image it does not do anything.

Am I missing something?

Thank you very much for your time,
Cristina

Why are you resampling? And if you are resampling, you should pay attention to spatial extent of your images, their overlap etc.

You could take a look at a tutorial from earlier this year, or its video recording. Of particular interest would be ITK image grid notebook.

Thank you very much!

I actually donā€™t know why I was resampling. I finally managed to fix the problem by changing the
ResampleImageFilter to ChangeInformationImageFilter.

1 Like