Skip to content

Commit

Permalink
BUG: Fix and test casting of resampled label map.
Browse files Browse the repository at this point in the history
Change to a vtkImageCast to simplify the setting of the resampled label type.
Cleaned up the logic, removed an extraneous AddNode (the scalar volume reader
adds the input node to the scene).
Removed debugging print outs.
Removed the try block as it was masking errors.
Added two tests for the converter script, with new label map volumes that were
created from the test DICOM series. The second one has a rotation applied to it
to trigger the reformatting.
Use sys.exit to ensure that tests will be detected as failing.

Issue QIICR#42
  • Loading branch information
Nicole Aucoin committed Mar 25, 2015
1 parent 6810cd7 commit 3b8c582
Show file tree
Hide file tree
Showing 6 changed files with 144 additions and 100 deletions.
Binary file added Prototype/TestData/labelToSEG.nrrd
Binary file not shown.
Binary file added Prototype/TestData/labelToSEGTransformed.nrrd
Binary file not shown.
52 changes: 43 additions & 9 deletions Testing/Python/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -3,20 +3,22 @@
# Install test data
#
message("Module source dir is ${${MODULE_NAME}_SOURCE_DIR}")
set(TEST_DATA_DIR ${${MODULE_NAME}_SOURCE_DIR}/Prototype/TestData/DICOM.CT)
set(TEMP "${CMAKE_BINARY_DIR}/Testing/Temporary/DICOM.CT")
MESSAGE(STATUS "Reporting: copying files from ${TEST_DATA_DIR} to ${TEMP}")
set(REPORTING_TEST_DATA_DIR "${${MODULE_NAME}_SOURCE_DIR}/Prototype/TestData")
set(SOURCE_DICOM_DIR ${REPORTING_TEST_DATA_DIR}/DICOM.CT)
set(TEMP_DIR "${CMAKE_BINARY_DIR}/Testing/Temporary")
set(TEMP_DICOM "${TEMP_DIR}/DICOM.CT")
MESSAGE(STATUS "Reporting: copying files from ${SOURCE_DICOM_DIR} to ${TEMP_DICOM}")
configure_file(
${TEST_DATA_DIR}/instance_487.dcm
${TEMP}/instance_487.dcm
${SOURCE_DICOM_DIR}/instance_487.dcm
${TEMP_DICOM}/instance_487.dcm
COPYONLY)
configure_file(
${TEST_DATA_DIR}/instance_488.dcm
${TEMP}/instance_488.dcm
${SOURCE_DICOM_DIR}/instance_488.dcm
${TEMP_DICOM}/instance_488.dcm
COPYONLY)
configure_file(
${TEST_DATA_DIR}/instance_489.dcm
${TEMP}/instance_489.dcm
${SOURCE_DICOM_DIR}/instance_489.dcm
${TEMP_DICOM}/instance_489.dcm
COPYONLY)


Expand All @@ -37,3 +39,35 @@ slicerMacroBuildScriptedModule(
SCRIPTS "ReportingSelfTest.py"
RESOURCES ""
)

# test the converter script with this DICOM as well as label map volumes
# made from the DICOM volume
set(REPORTING_OUTPUT_SEG_DIR "${TEMP_DIR}/OutputSeg")
message("Using output segmentation dir of ${REPORTING_OUTPUT_SEG_DIR}")
file(MAKE_DIRECTORY ${REPORTING_OUTPUT_SEG_DIR})
message("Copying label volumes from ${REPORTING_TEST_DATA_DIR} to ${TEMP_DIR}")
configure_file(
${REPORTING_TEST_DATA_DIR}/labelToSEG.nrrd
${TEMP_DIR}/labelToSEG.nrrd
COPYONLY)
configure_file(
${REPORTING_TEST_DATA_DIR}/labelToSEGTransformed.nrrd
${TEMP_DIR}/labelToSEGTransformed.nrrd
COPYONLY)

add_test(
NAME "py_LabelToDICOMSEGConverterScript_TestLabelVolume"
COMMAND ${Slicer_LAUNCHER_EXECUTABLE} --no-main-window
--additional-module-path ${CMAKE_BINARY_DIR}/${Slicer_QTLOADABLEMODULES_LIB_DIR}
--python-script ${${MODULE_NAME}_SOURCE_DIR}/Util/LabelToDICOMSEGConverterScript.py
${TEMP_DICOM} ${TEMP_DIR}/labelToSEG.nrrd ${REPORTING_OUTPUT_SEG_DIR}
)
# add a test that triggers the resampling of the input label map
add_test(
NAME "py_LabelToDICOMSEGConverterScript_TestTransformedLabelVolume"
COMMAND ${Slicer_LAUNCHER_EXECUTABLE} --no-main-window
--additional-module-path ${CMAKE_BINARY_DIR}/${Slicer_QTLOADABLEMODULES_LIB_DIR}
--python-script ${${MODULE_NAME}_SOURCE_DIR}/Util/LabelToDICOMSEGConverterScript.py
${TEMP_DICOM} ${TEMP_DIR}/labelToSEGTransformed.nrrd ${REPORTING_OUTPUT_SEG_DIR}
)

2 changes: 1 addition & 1 deletion Testing/Python/ReportingSelfTest.py
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,7 @@ def __init__(self, parent):
except AttributeError:
slicer.selfTests = {}
slicer.selfTests['ReportingSelfTest'] = self.runTest
print slicer.selfTests
# print 'ReportingSelfTest: Full list of Slicer self tests = ',slicer.selfTests

def runTest(self):
tester = ReportingSelfTestTest()
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,7 @@ def __init__(self, parent):
except AttributeError:
slicer.selfTests = {}
slicer.selfTests['LabelToDICOMSEGConverterSelfTest'] = self.runTest
print slicer.selfTests
# print 'LabelToDICOMSEGConverterSelfTest: full list of Slicer self tests = ',slicer.selfTests

def runTest(self):
tester = LabelToDICOMSEGConverterSelfTestTest()
Expand Down
188 changes: 99 additions & 89 deletions Util/LabelToDICOMSEGConverterScript.py
Original file line number Diff line number Diff line change
Expand Up @@ -51,6 +51,11 @@
def DoIt(inputDir, labelFile, outputDir, forceLabel):

dbDir1 = slicer.app.temporaryPath+'/LabelConverter'

if not hasattr(slicer.modules, 'reporting'):
print 'The Reporting module has not been loaded into Slicer, script cannot run!\n\tTry setting the --additional-module-path parameter.'
sys.exit(1)

reportingLogic = slicer.modules.reporting.logic()

print('Temporary directory location: '+dbDir1)
Expand All @@ -60,52 +65,50 @@ def DoIt(inputDir, labelFile, outputDir, forceLabel):
if slicer.dicomDatabase:
dbDir0 = os.path.split(slicer.dicomDatabase.databaseFilename)[0]

try:
dicomWidget = slicer.modules.dicom.widgetRepresentation().self()
dicomWidget.onDatabaseDirectoryChanged(dbDir1)

# import DICOM study
indexer = ctk.ctkDICOMIndexer()
indexer.addDirectory(slicer.dicomDatabase, inputDir, None)
indexer.waitForImportFinished()

print('DICOM import finished!')

#
# Read the input DICOM series as a volume
#
dcmList = []
for dcm in os.listdir(inputDir):
if len(dcm)-dcm.rfind('.dcm') == 4:
dcmList.append(inputDir+'/'+dcm)

scalarVolumePlugin = slicer.modules.dicomPlugins['DICOMScalarVolumePlugin']()

loadables = scalarVolumePlugin.examine([dcmList])

if len(loadables) == 0:
print 'Could not parse the DICOM Study!'
exit()

inputVolume = scalarVolumePlugin.load(loadables[0])
slicer.mrmlScene.AddNode(inputVolume)
print('Input volume loaded!')

# read the label volume
labelVolume = slicer.vtkMRMLScalarVolumeNode()
sNode = slicer.vtkMRMLVolumeArchetypeStorageNode()
sNode.SetFileName(labelFile)
sNode.ReadData(labelVolume)
labelVolume.LabelMapOn()

if forceLabel>0:
print('Forcing label to '+str(forceLabel))
labelImage = labelVolume.GetImageData()
thresh = vtk.vtkImageThreshold()
if vtk.vtkVersion().GetVTKMajorVersion() < 6:
thresh.SetInput(labelImage)
else:
thresh.SetInputData(labelImage)
dicomWidget = slicer.modules.dicom.widgetRepresentation().self()
dicomWidget.onDatabaseDirectoryChanged(dbDir1)

# import DICOM study
indexer = ctk.ctkDICOMIndexer()
indexer.addDirectory(slicer.dicomDatabase, inputDir, None)
indexer.waitForImportFinished()

print('DICOM import finished!')

#
# Read the input DICOM series as a volume
#
dcmList = []
for dcm in os.listdir(inputDir):
if len(dcm)-dcm.rfind('.dcm') == 4:
dcmList.append(inputDir+'/'+dcm)

scalarVolumePlugin = slicer.modules.dicomPlugins['DICOMScalarVolumePlugin']()

loadables = scalarVolumePlugin.examine([dcmList])

if len(loadables) == 0:
print 'Could not parse the DICOM Study!'
sys.exit(1)

inputVolume = scalarVolumePlugin.load(loadables[0])
print 'Input volume loaded! ID = ', inputVolume.GetID()

# read the label volume
labelVolume = slicer.vtkMRMLScalarVolumeNode()
sNode = slicer.vtkMRMLVolumeArchetypeStorageNode()
sNode.SetFileName(labelFile)
sNode.ReadData(labelVolume)
labelVolume.LabelMapOn()

if forceLabel>0:
# print('Forcing label to '+str(forceLabel))
labelImage = labelVolume.GetImageData()
thresh = vtk.vtkImageThreshold()
if vtk.vtkVersion().GetVTKMajorVersion() < 6:
thresh.SetInput(labelImage)
else:
thresh.SetInputData(labelImage)
thresh.ThresholdBetween(1, labelImage.GetScalarRange()[1])
thresh.SetInValue(int(forceLabel))
thresh.SetOutValue(0)
Expand All @@ -115,47 +118,54 @@ def DoIt(inputDir, labelFile, outputDir, forceLabel):
labelImage = thresh.GetOutput()
labelVolume.SetAndObserveImageData(labelImage)

slicer.mrmlScene.AddNode(labelVolume)

volumesLogic = slicer.modules.volumes.logic()
geometryCheckString = volumesLogic.CheckForLabelVolumeValidity(inputVolume, labelVolume)
if geometryCheckString != "":
print('Label volume geometry mismatch, resampling:\n%s' % geometryCheckString)

# resample label to the input volume raster
resampledLabel = slicer.vtkMRMLScalarVolumeNode()
slicer.mrmlScene.AddNode(resampledLabel)
resampledLabel = volumesLogic.ResampleVolumeToReferenceVolume(labelVolume, inputVolume)
if resampledLabel.GetImageData() != None:
# double check that the pixel type is unsigned short
if resampledLabel.GetImageData().GetScalarType() != vtk.VTK_UNSIGNED_SHORT:
if vtk.vtkVersion().GetVTKMajorVersion() < 6:
resampledLabel.GetImageData().SetScalarType(vtk.VTK_UNSIGNED_SHORT)
else:
tp = vtk.vtkTrivialProducer()
tp.SetOutput(resampledLabel.GetImageData())
outInfo = tp.GetOutputInformation(0)
vtk.vtkDataObject().SetPointDataActiveScalarInfo(outInfo, vtk.VTK_UNSIGNED_SHORT, vtk.vtkImageData().GetNumberOfScalarComponents(outInfo))
labelVolume = resampledLabel

displayNode = slicer.vtkMRMLLabelMapVolumeDisplayNode()
displayNode.SetAndObserveColorNodeID(reportingLogic.GetDefaultColorNode().GetID())
slicer.mrmlScene.AddNode(displayNode)

labelVolume.SetAttribute('AssociatedNodeID',inputVolume.GetID())
labelVolume.LabelMapOn()
labelVolume.SetAndObserveDisplayNodeID(displayNode.GetID())

# initialize the DICOM DB for Reporting logic, save as DICOM SEG
labelCollection = vtk.vtkCollection()
labelCollection.AddItem(labelVolume)

print('About to write DICOM SEG!')
dbFileName = slicer.dicomDatabase.databaseFilename
reportingLogic.InitializeDICOMDatabase(dbFileName)
reportingLogic.DicomSegWrite(labelCollection, outputDir)
except:
print('Error occurred!')
slicer.mrmlScene.AddNode(labelVolume)
print 'Label volume added, id = ', labelVolume.GetID()

volumesLogic = slicer.modules.volumes.logic()
geometryCheckString = volumesLogic.CheckForLabelVolumeValidity(inputVolume, labelVolume)
if geometryCheckString != "":
# resample label to the input volume raster
resampledLabel = slicer.vtkMRMLScalarVolumeNode()
slicer.mrmlScene.AddNode(resampledLabel)
print 'Resampled label added, id = ', resampledLabel.GetID()
resampledLabel = volumesLogic.ResampleVolumeToReferenceVolume(labelVolume, inputVolume)
if resampledLabel.GetImageData() != None:
# double check that the pixel type is unsigned short
if resampledLabel.GetImageData().GetScalarType() != vtk.VTK_UNSIGNED_SHORT:
cast = vtk.vtkImageCast()
cast.SetOutputScalarTypeToUnsignedShort()
# set the label volume to be the cast resampled label volume
if vtk.vtkVersion().GetVTKMajorVersion() < 6:
cast.SetInput(resampledLabel.GetImageData())
cast.Update()
labelVolume.SetAndObserveImageData(cast.GetOutput())
else:
cast.SetInputConnection(resampledLabel.GetImageDataConnection())
cast.Update()
labelVolume.SetImageDataConnection(cast.GetOutputPort())
if labelVolume.GetImageData().GetScalarType() != vtk.VTK_UNSIGNED_SHORT:
print 'Failed to cast label volume to unsigned short, type is ', vtk.vtkImageScalarTypeNameMacro(labelVolume.GetImageData().GetScalarType())
sys.exit(1)
else:
# the resampled label map is the right type, use it
labelVolume = resampledLabel

displayNode = slicer.vtkMRMLLabelMapVolumeDisplayNode()
displayNode.SetAndObserveColorNodeID(reportingLogic.GetDefaultColorNode().GetID())
slicer.mrmlScene.AddNode(displayNode)

labelVolume.SetAttribute('AssociatedNodeID',inputVolume.GetID())
labelVolume.LabelMapOn()
labelVolume.SetAndObserveDisplayNodeID(displayNode.GetID())

# initialize the DICOM DB for Reporting logic, save as DICOM SEG
labelCollection = vtk.vtkCollection()
labelCollection.AddItem(labelVolume)

print('About to write DICOM SEG!')
dbFileName = slicer.dicomDatabase.databaseFilename
reportingLogic.InitializeDICOMDatabase(dbFileName)
reportingLogic.DicomSegWrite(labelCollection, outputDir)

dicomWidget.onDatabaseDirectoryChanged(dbDir0)

Expand All @@ -164,12 +174,12 @@ def DoIt(inputDir, labelFile, outputDir, forceLabel):
if len(sys.argv)<4:
print 'Input parameters missing!'
print 'Usage: ',sys.argv[0],' <input directory with DICOM data> <input label> <output dir> <optional: label to assign>'
exit()
sys.exit(1)
else:
inputDICOMDir = sys.argv[1]
inputLabelName = sys.argv[2]
outputDICOMDir = sys.argv[3]
forceLabel = 0
if len(sys.argv)>4:
forceLabel = sys.argv[4]
DoIt(inputDICOMDir, inputLabelName, outputDICOMDir, forceLabel)
DoIt(inputDICOMDir, inputLabelName, outputDICOMDir, forceLabel)

0 comments on commit 3b8c582

Please sign in to comment.