MHA Numpy

I’m trying to write this code in Python, but i don’t know what i have to do when i have to create an iterator to build the mask, because i don’t know what are the values of the array created from mha file.

typedef itk::ImageRegionIterator IteratorType;
IteratorType maskIt(outputImage,outputImage->GetLargestPossibleRegion());
LabelPixelType backgroundValue = 0;
LabelPixelType foregroundValue = 1;
for (maskIt.GoToBegin(); !maskIt.IsAtEnd(); ++maskIt)
{
ImageType::PointType point;
outputImage->TransformIndexToPhysicalPoint(maskIt.GetIndex(),point);
if (maskInterpolator->IsInsideBuffer(point))
{
if (maskInterpolator->Evaluate(point) == 0)
{
maskIt.Set(backgroundValue);
}
else
{
maskIt.Set(foregroundValue);
}
}
else
{
maskIt.Set(backgroundValue);
}
}

Image iterator are not available from Python, because accessing images via this mechanism is too slow due to wrapping overhead.

I’m using NumPy to access the file, but i don’t understand what represents elements in the array that i get with array = itk.GetArrayViewFromImage(outputImage)

Array elements are intensities of pixels in the image. Try reading a super small image (lets say 3x4x5) with distinctive values, and then access that via itk.GetArrayViewFromImage.

1 Like

Thank you so much, but when i have to use TransformIndexToPhysicalPoint what i have to put to replace maskIt.GetIndex() ?
EDIT: the first value of the matrix is the slice number, the second is the y axe, the third is the x axe and all the elements cointaned in the matrix are the pixel values

Why don’t you just resample the Mask image onto the output image with the ResampleImageFilter, then apply a BinaryThresholdImageFilter?

3 Likes

Update:

pixelIndex=itk.Index[2] ()
for z in range (array.shape[0]):#slice
for y in range (array.shape[1]): #Y axes
for x in range (array.shape[2]): #X axes
pixelIndex[0]=x
pixelIndex[1]=y
point=ImageType.New()
outputImage.TransformIndexToPhysicalPoint(pixelIndex,point)
if (maskInterpolator.IsInsideBuffer(point)):
if (maskInterpolator.Evaluate(point) == 0):
maskIt.Set(backgroundValue)
else:
maskIt.Set(foregroundValue)
else:
maskIt.Set(backgroundValue)

I obtain this error:

TypeError: itkImageBase3_TransformIndexToPhysicalPoint expected 2 arguments, got 3

This may be helpful:

Dim=3
PixelType=itk.F
LabelPixelType=itk. US
ImageType=itk.Image[PixelType,Dim]

@blowekamp and @dzenanz have good suggestions – the result will be simpler and faster if itk.resample_image_filter with use_reference_image=True and reference_image=input_image followed by a call to itk.binary_threshold_image_filter.

See also:

1 Like

Thank you so much, i’ll follow your tips, but i would like to understand also the previous error. Do you have any idea?

It’s a little bit different from your interpretation because i have to check if a pixel is inside or outside a ROI contained in a TIF as an image sequence. I have already created an interpolation with itk.NearestNeighborInterpolateImageFunction so now i have to check if a point is inside or outside to create the mask

There are two options for TransformIndexToPhysicalPoint, one that accepts a reference to the Point for its output, and one that returns the Point as its output. References are not supported well in Python, so you want the one that returns the value. Therefore,

point = outputImage.TransformIndexToPhysicalPoint(pixelIndex)

A call to itk.mask_image_filter may also be helpful:

1 Like

I followed your suggestions but i have a problem because the mask and the input image don’t occupy the same physical space, despite having assigned them identical size,origin etc. It seems that the value of the mask resets in the function MaskImageFilter

maskImage=flipFilter.GetOutput()
maskImage.SetOrigin(reader4mm.GetOutput().GetOrigin()) >>>Assignation values
maskImage.SetSpacing(reader4mm.GetOutput().GetSpacing())
maskImage.SetDirection(reader4mm.GetOutput().GetDirection())
#create the final mask
MaskFilterType=itk.MaskImageFilter[LabeledImageType,MaskImageType,LabeledImageType]
maskFilter=MaskFilterType.New()
maskFilter.SetInput(outputImage)
maskFilter.SetMaskImage(maskImage) >>>Mask seems to reset values
WriterType=itk.ImageFileWriter[LabeledImageType]
writer=WriterType.New()
writer.SetFileName(path+"/mask.mha")
writer.SetInput(maskFilter.GetOutput())
try:
writer.Update();
except Exception as error:
print(“EXIT FAILURE”)
print(error)
sys.exit(1)
print(“EXIT SUCCESS”)

When the pipeline runs, the flipFilter can reset the value of the output – use a ChangeInformationImageFilter in this case:

https://itk.org/Doxygen/html/classitk_1_1ChangeInformationImageFilter.html

It’s possible that ChangeInformationImageFilter is not implemented in python? I tried to do something like this:

maskImage=flipFilter.GetOutput()
origin=reader4mm.GetOutput().GetOrigin()
spacing=reader4mm.GetOutput().GetSpacing()
direction=reader4mm.GetOutput().GetDirection()
maskImage.SetOutputSpacing(spacing)
maskImage.ChangeSpacingOn()

but it says that AttributeError: 'itkImageUC3' object has no attribute 'ChangeSpacingOn'.
It may be helpful to say that i used flipFilter (itk.FlipImageFilter) to change the axis of the mask and that the mask should have the same size of the output image [256,256,50] (so i have to set: origin, direction and spacing)

EDIT: i tried to not flip the axes with flipFilter and to mantain

maskImage.SetOrigin(reader4mm.GetOutput().GetOrigin())
maskImage.SetSpacing(reader4mm.GetOutput().GetSpacing())
maskImage.SetDirection(reader4mm.GetOutput().GetDirection())

but when i use the function MaskImageFilter the values of maskImage (origin,direction,spacing) seems to be reset in the function:

c:\p\ipp\standalone-build\itks\modules\core\common\include\itkImageToImageFilter.hxx:220:
itk::ERROR: MaskImageFilter(000002A7B7568450): Inputs do not occupy the same physical space!
InputImage Origin: [-1.6488900e+02, -7.2781100e+01, 1.9833500e+02], InputImage_1 Origin: [0.0000000e+00, 0.0000000e+00, 4.9000000e+01]
Tolerance: 1.4063000e-06
InputImage Spacing: [1.4063000e+00, 1.4063000e+00, 4.0000020e+00], InputImage_1 Spacing: [1.0000000e+00, 1.0000000e+00, 1.0000000e+00]
Tolerance: 1.4063000e-06
InputImage Direction: 1.0000000e+00 0.0000000e+00 0.0000000e+00
0.0000000e+00 -0.0000000e+00 1.0000000e+00
0.0000000e+00 -1.0000000e+00 -0.0000000e+00
, InputImage_1 Direction: 1.0000000e+00 0.0000000e+00 0.0000000e+00
0.0000000e+00 1.0000000e+00 0.0000000e+00
0.0000000e+00 0.0000000e+00 -1.0000000e+00

but if i print maskImage in the code the values of origin, dimension and direction are corrected

ChangeInformationImageFilter is wrapped. It could be called like this:

read_image = reader4mm.GetOutput()
maskImage = itk.change_information_image_filter(flipFilter.GetOutput(),
    output_origin=read_image.GetOrigin(),
    output_spacing=read_image.GetSpacing(),
    output_direction=read_image.GetDirection()
)

But, depending on what you are trying to do and your data, you may want / need to use ResampleImageFilter as suggested earlier.

1 Like