I have a small python VTK
function that calculates the volume and surface area of an object embedded in a stack of TIFF
images. To read the TIFF's
into VTK
, I have used vtkTIFFReader
and processed the result using vtkImageThreshold
. I then use vtkMassProperties
to extract the volume and surface area of the object identified after thresholding.
With VTK-5.04
, this function returns the correct value for a test stack (3902 pixels). However, using VTK-5.4.2
the same function returns a different value (422 pixels). Can someone explain this?
def testvtk():
# read 36 TIFF images. Each TIFF is 27x27 pixels
v16=vtk.vtkTIFFReader()
v16.SetFilePrefix("d:/test/slice")
v16.SetDataExtent(0,27,0,27,1,36)
v16.SetFilePattern("%s%04d.tif")
v16.SetDataSpacing (1,1,1)
v16.Update()
# Threshold level for seperating background/foreground pixels
maxthres=81
# Threshold the image stack
thres=vtk.vtkImageThreshold()
thres.SetInputConnection(v16.GetOutputPort())
thres.ThresholdByLower(0)
thres.ThresholdByUpper(maxthres)
# create ISO surface from thresholded images
iso=vtk.vtkImageMarchingCubes()
iso.SetInputConnection(thres.GetOutputPort())
# Have VTK calculate the Mass (volume) and surface area
Mass = vtk.vtkMassProperties()
Mass.SetInputConnection(iso.GetOutputPort())
Mass.Update()
# just print the results
print "Volume = ", Mass.GetVolume()
print "Surface = ", Mass.GetSurfaceArea()
By testing both VTK-5.4.2 and VTK-5.2.1, I narrowed things down a bit and believe this behaviour was introduced between versions 5.0.4 and 5.2.1.
It seems that in VTK-5.4.2, vtkTIFFReader ignores the x and y values set in the SetDataSpacing method. Instead vtkTIFFReader is calculating the x and y dataspacing from the resolution reported by the TIFF files.
I have never heard of VTK before, but here it goes.
Good thing about opensource software is that you can check the source code directly. Better yet, if there's web-based version control browser, we can talk about it online like this.
Let's see vtkMassProperties
in question. 5.0.4 uses r1.28 and 5.4.2 uses r1.30. Here's the diff between r1.28 and r.30. The part that may affect volume calculations are
vol[2] += (area * (double)u[2] * (double)zavg); // 5.0.4
vol[2] += (area * u[2] * zavg); // 5.4.2
and
kxyz[0] = (munc[0] + (wxyz/3.0) + ((wxy+wxz)/2.0)) /(double)(numCells); // 5.0.4
kxyz[0] = (munc[0] + (wxyz/3.0) + ((wxy+wxz)/2.0)) /numCells; // 5.4.2
but all changes look ok to me.
Next suspicious are the vtkMarchingCubes
. Diff between r1.1.6.1 and 1.5.
self->UpdateProgress ((double) k / ((double) dims[2] - 1)); // 5.0.4
self->UpdateProgress (k / static_cast<double>(dims[2] - 1)); // 5.4.2
and
estimatedSize = (int) pow ((double) (dims[0] * dims[1] * dims[2]), .75); // 5.0.4
estimatedSize = static_cast<int>(
pow(static_cast<double>(dims[0]*dims[1]*dims[2]),0.75)); // 5.4.2
Again, they are fixing stuff on casting, but it looks ok.
Might as well see vtkImageThreshold
too. Diff between r1.50 and r1.52.
lowerThreshold = (IT) inData->GetScalarTypeMin(); // 5.0.4
lowerThreshold = static_cast<IT>(inData->GetScalarTypeMin()); // 5.4.2
There are bunch more, but they are all casting stuff.
It gets more interesting with vtkTIFFReader
. Diff between 1.51 and 1.63. As you can see by the difference in the revision numbers, there has been some development in this class compared to others. Here are the checkin comments:
From the amount of changes that was made in vtkTIFFReader, I would guess that the difference in behavior is coming from there. For example, it may have started to recognize your Tiff as different format and changed the internal pixel values. Try printing out pixel values and see if there is any difference. If the pixel values have changed maxthres=81
may be too high.
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With