Negative Values in DVF Histogram

I tried to get the histogram for each segment with this code.
However, I got the same result for every segment.
I apologize for the inconvenience. Please check.

>>> labelValue = 1
Loading with imageIOName: GDCM
>>> labelmapVolumeNode =slicer.mrmlScene.AddNewNodeByClass("vtkMRMLLabelMapVolumeNode")
>>> segmentationNode = getNode('segmentation')
>>> volumeNode = getNode('Displacement Field')
>>> labelmapVolumeNode = getNode('LabelMapVolume')
>>> slicer.modules.segmentations.logic().ExportVisibleSegmentsToLabelmapNode(segmentationNode, labelmapVolumeNode, volumeNode)
True
>>> volumeArray = slicer.util.arrayFromVolume(volumeNode)
>>> labelArray = slicer.util.arrayFromVolume(labelmapVolumeNode)
>>> segmentVoxels = volumeArray[labelArray==labelValue]
>>> import numpy as np
>>> histogram = np.histogram(arrayFromVolume(volumeNode), bins=100)
>>> slicer.util.plot(histogram, xColumnIndex = 1)
(MRMLCorePython.vtkMRMLPlotChartNode)000002174751D408
>>> range = [-5.0, 5.0]
>>> bins = 100
>>> histogramR, histogramBins = np.histogram(arrayFromVolume(volumeNode)[:,:,:,0], bins, range)
>>> histogramA = np.histogram(arrayFromVolume(volumeNode)[:,:,:,1], bins, range)[0]
>>> histogramS = np.histogram(arrayFromVolume(volumeNode)[:,:,:,2], bins, range)[0]
>>> slicer.util.plot(np.vstack([histogramBins[:-1], histogramR, histogramA, histogramS]).T, xColumnIndex = 0, columnNames=['N', 'R', 'A', 'S'], title='Histogram')
(MRMLCorePython.vtkMRMLPlotChartNode)000002174751A3A8
>>> magnitudes = np.linalg.norm(arrayFromVolume(volumeNode), axis=3)
>>> histogram = np.histogram(magnitudes, bins=100)
>>> slicer.util.plot(histogram, xColumnIndex = 1, columnNames=['magnitude', 'N'], title="displacement")
(MRMLCorePython.vtkMRMLPlotChartNode)000002174751AE88
>>>

You nicely computed segmentVoxels but did not use it when you computed the histogram.

segmentVoxels = volumeArray[labelArray==labelValue]

Is the above not correct?
If you know, I would like to know more about the code to add and where to fix it.

labelValue = 1
Loading with imageIOName: GDCM
labelmapVolumeNode =slicer.mrmlScene.AddNewNodeByClass(ā€œvtkMRMLLabelMapVolumeNodeā€)
segmentationNode = getNode(ā€˜segmentationā€™)
volumeNode = getNode(ā€˜Displacement Fieldā€™)
labelmapVolumeNode = getNode(ā€˜LabelMapVolumeā€™)
slicer.modules.segmentations.logic().ExportVisibleSegmentsToLabelmapNode(segmentationNode, labelmapVolumeNode, volumeNode)
True
volumeArray = slicer.util.arrayFromVolume(volumeNode)
labelArray = slicer.util.arrayFromVolume(labelmapVolumeNode)
segmentVoxels = volumeArray[labelArray==labelValue]
import numpy as np
histogram = np.histogram(arrayFromVolume(volumeNode),segmentVoxels, bins=50)
Traceback (most recent call last):
File ā€œā€, line 1, in
File ā€œ<array_function internals>ā€, line 4, in histogram
TypeError: _histogram_dispatcher() got multiple values for argument ā€˜binsā€™

I changed the code a bit.

histogram = np.histogram(arrayFromVolume(volumeNode), bins=100)

Added a segment voxel to the above.

histogram = np.histogram(arrayFromVolume(volumeNode),segmentVoxels, bins=100)

But an error occurred.
Iā€™m sorry for the inconvenience, but please confirm.

Sorry for asking so many questions.
Next, I tried to write such a code.
I changed the arrayFromVolume(volumeNode) in the code you gave me to segmentVoxels.
However, I got an error again.
How can I get the histogram for each segment?

>>> labelValue = 1
Loading with imageIOName: GDCM
>>> labelmapVolumeNode =slicer.mrmlScene.AddNewNodeByClass("vtkMRMLLabelMapVolumeNode")
>>> segmentationNode = getNode('segmentation')
>>> volumeNode = getNode('Displacement Field')
>>> labelmapVolumeNode = getNode('LabelMapVolume')
>>> slicer.modules.segmentations.logic().ExportVisibleSegmentsToLabelmapNode(segmentationNode, labelmapVolumeNode, volumeNode)
True
>>> volumeArray = slicer.util.arrayFromVolume(volumeNode)
>>> labelArray = slicer.util.arrayFromVolume(labelmapVolumeNode)
>>> import numpy as np
>>> histogram = np.histogram(arrayFromVolume(volumeNode), bins=100)
>>> slicer.util.plot(histogram, xColumnIndex = 1)
(MRMLCorePython.vtkMRMLPlotChartNode)000001F95ACA9E88
>>> segmentationNode = getNode('segmentation')
>>> volumeNode = getNode('Displacement Field')
>>> labelmapVolumeNode = getNode('LabelMapVolume')
>>> slicer.modules.segmentations.logic().ExportVisibleSegmentsToLabelmapNode(segmentationNode, labelmapVolumeNode, volumeNode)
True
>>> volumeArray = slicer.util.arrayFromVolume(volumeNode)
>>> labelArray = slicer.util.arrayFromVolume(labelmapVolumeNode)
>>> segmentVoxels = volumeArray[labelArray==labelValue]
>>> import numpy as np
>>> histogram = np.histogram(arrayFromVolume(volumeNode), bins=100)
>>> slicer.util.plot(histogram, xColumnIndex = 1)
(MRMLCorePython.vtkMRMLPlotChartNode)000001F95ACA9408
>>> segmentationNode = getNode('segmentation')
>>> volumeNode = getNode('Displacement Field')
>>> labelmapVolumeNode = getNode('LabelMapVolume')
>>> slicer.modules.segmentations.logic().ExportVisibleSegmentsToLabelmapNode(segmentationNode, labelmapVolumeNode, volumeNode)
True
>>> volumeArray = slicer.util.arrayFromVolume(volumeNode)
>>> labelArray = slicer.util.arrayFromVolume(labelmapVolumeNode)
>>> segmentVoxels = volumeArray[labelArray==labelValue]
>>> import numpy as np
>>> histogram = np.histogram(segmentVoxels, bins=100)
>>> slicer.util.plot(histogram, xColumnIndex = 1)
(MRMLCorePython.vtkMRMLPlotChartNode)000001F95ACA9E28
>>> range = [-5.0, 5.0]
>>> bins = 100
>>> histogramR, histogramBins = np.histogram(segmentVoxles[:,:,:,0], bins, range)
Traceback (most recent call last):
File "<console>", line 1, in <module>
NameError: name 'segmentVoxles' is not defined

There is a typo in the variable name.

Where is the mistake?
I changed arrayFromVolume(volumeNode) to segmentVoxels, but this did not give me a histogram for each segment.
What should I do?

labelValue = 1
Loading with imageIOName: GDCM
labelmapVolumeNode =slicer.mrmlScene.AddNewNodeByClass(ā€œvtkMRMLLabelMapVolumeNodeā€)
segmentationNode = getNode(ā€˜segmentationā€™)
volumeNode = getNode(ā€˜Displacement Fieldā€™)
labelmapVolumeNode = getNode(ā€˜LabelMapVolumeā€™)
slicer.modules.segmentations.logic().ExportVisibleSegmentsToLabelmapNode(segmentationNode, labelmapVolumeNode, volumeNode)
True
volumeArray = slicer.util.arrayFromVolume(volumeNode)
labelArray = slicer.util.arrayFromVolume(labelmapVolumeNode)
segmentVoxels = volumeArray[labelArray==labelValue]
import numpy as np
histogram = np.histogram(segmentVoxels, bins=100)
slicer.util.plot(histogram, xColumnIndex = 1)
(MRMLCorePython.vtkMRMLPlotChartNode)0000024502908828
histogram = np.histogram(arrayFromVolume(volumeNode), bins=100)
slicer.util.plot(histogram, xColumnIndex = 1)
(MRMLCorePython.vtkMRMLPlotChartNode)0000024502908408
range = [-5.0, 5.0]
bins = 100
histogramR, histogramBins = np.histogram(segmentationVoxels[:,:,:,0], bins, range)
Traceback (most recent call last):
File ā€œā€, line 1, in
NameError: name ā€˜segmentationVoxelsā€™ is not defined

If you know how to get each contour, please let me know.

segmentVoxles is mistyped. The name is not the same as the segmentVoxels variable that you created above.

Thank you for letting me know.
Iā€™ve modified the code.
Is it now possible to get the DVF histogram for each segment?

labelValue = 1
Loading with imageIOName: GDCM
labelmapVolumeNode =slicer.mrmlScene.AddNewNodeByClass(ā€œvtkMRMLLabelMapVolumeNodeā€)
segmentationNode = getNode(ā€˜segmentationā€™)
volumeNode = getNode(ā€˜Displacement Fieldā€™)
labelmapVolumeNode = getNode(ā€˜LabelMapVolumeā€™)
slicer.modules.segmentations.logic().ExportVisibleSegmentsToLabelmapNode(segmentationNode, labelmapVolumeNode, volumeNode)
True
volumeArray = slicer.util.arrayFromVolume(volumeNode)
labelArray = slicer.util.arrayFromVolume(labelmapVolumeNode)
segmentVoxels = volumeArray[labelArray==labelValue]
import numpy as np
histogram = np.histogram(segmentVoxels, bins=100)
slicer.util.plot(histogram, xColumnIndex = 1)
(MRMLCorePython.vtkMRMLPlotChartNode)0000015B0A709EE8
range = [-5.0, 5.0]
bins = 100
histogramR, histogramBins = np.histogram(segmentVoxels[:,0], bins, range)
histogramA = np.histogram(segmentVoxels[:,1], bins, range)[0]
histogramS = np.histogram(segmentVoxels[:,2], bins, range)[0]
slicer.util.plot(np.vstack([histogramBins[:-1], histogramR, histogramA, histogramS]).T, xColumnIndex = 0, columnNames=[ā€˜Nā€™, ā€˜Rā€™, ā€˜Aā€™, ā€˜Sā€™], title=ā€˜Histogramā€™)
(MRMLCorePython.vtkMRMLPlotChartNode)0000015B0A709468
magnitudes = np.linalg.norm(segmentVoxels, axis=1)
histogram = np.histogram(magnitudes, bins=100)
slicer.util.plot(histogram, xColumnIndex = 1, columnNames=[ā€˜magnitudeā€™, ā€˜Nā€™], title=ā€œdisplacementā€)
(MRMLCorePython.vtkMRMLPlotChartNode)0000015B0DBC3F48

The code you gave me before gave me this error.

histogramR, histogramBins = np.histogram(segmentVoxels[:,:,:,0], bins, range)
Traceback (most recent call last):
File ā€œā€, line 1, in
IndexError: too many indices for array: array is 2-dimensional, but 4 were indexed

magnitudes = np.linalg.norm(segmentVoxels, axis=3)
Traceback (most recent call last):
File ā€œā€, line 1, in
File ā€œ<array_function internals>ā€, line 6, in norm
return sqrt(add.reduce(s, axis=axis, keepdims=keepdims))
numpy.AxisError: axis 3 is out of bounds for array of dimension 2

I modified this part, but is the code at the top correct?

The code snippets that I provided above all work correctly and you only need a little Python and numpy knowledge to apply them to your data. You may try to complete Python and numpy tutorials or find someone locally who can help you out with this.

Is it wrong to change [:,:,:,0] and axis=3 to fit the segment?
I was able to get the DVF for each axis and the total DVF with this code, but was the code wrong?
Is it impossible to get the DVF for each segment?
Is it possible to calculate only by volume?
I would like your answer.
image
image

You can use standard Python array indexing to filter voxels that are in a selected segment. For example:

volumeNode = getNode('Displacement Field')
segmentationNode = getNode('Segmentation')
segmentId = 'Segment_3'

magnitudes = np.linalg.norm(arrayFromVolume(volumeNode), axis=3)
magnitudesInSegment = magnitudes[arrayFromSegmentBinaryLabelmap(segmentationNode, segmentId) != 0]
histogram = np.histogram(magnitudesInSegment, bins=100)
slicer.util.plot(histogram, xColumnIndex = 1, columnNames=['magnitude', 'N'], title="displacement")

Thank you for your answer.
I tried, but I got an error.
I apologize for the inconvenience, but I would appreciate your advice.

>>> labelValue = 1
Loading with imageIOName: GDCM
>>> labelmapVolumeNode =slicer.mrmlScene.AddNewNodeByClass("vtkMRMLLabelMapVolumeNode")
>>> volumeNode = getNode('Displacement Field')
>>> segmentationNode = getNode('segmentation')
>>> segmentId = 'patient'
>>> labelmapVolumeNode = getNode('LabelMapVolume')
>>> slicer.modules.segmentations.logic().ExportVisibleSegmentsToLabelmapNode(segmentationNode, labelmapVolumeNode, volumeNode)
True
>>> volumeArray = slicer.util.arrayFromVolume(volumeNode)
>>> labelArray = slicer.util.arrayFromVolume(labelmapVolumeNode)
>>> segmentVoxels = volumeArray[labelArray==labelValue]
>>> import numpy as np
>>> histogram = np.histogram(arrayFromVolume(volumeNode), bins=100)
>>> slicer.util.plot(histogram, xColumnIndex = 1)
(MRMLCorePython.vtkMRMLPlotChartNode)00000222A95B87C8
>>> range = [-5.0, 5.0]
>>> bins = 100
>>> histogramR, histogramBins = np.histogram(arrayFromVolume(volumeNode)[:,:,:,0], bins, range)
>>> histogramA = np.histogram(arrayFromVolume(volumeNode)[:,:,:,1], bins, range)[0]
>>> histogramS = np.histogram(arrayFromVolume(volumeNode)[:,:,:,2], bins, range)[0]
>>> slicer.util.plot(np.vstack([histogramBins[:-1], histogramR, histogramA, histogramS]).T, xColumnIndex = 0, columnNames=['N', 'R', 'A', 'S'], title='Histogram')
(MRMLCorePython.vtkMRMLPlotChartNode)00000222A95B86A8
>>> magnitudes = np.linalg.norm(arrayFromVolume(volumeNode), axis=3)
>>> magnitudesInSegment = magnitudes[arrayFromSegmentBinaryLabelmap(segmentationNode, segmentId) != 0]
Traceback (most recent call last):
File "<console>", line 1, in <module>
IndexError: boolean index did not match indexed array along dimension 0; dimension is 128 but corresponding boolean dimension is 130

Iā€™ve tried drawing other code, but I still get this error.
I apologize for the inconvenience, but I would appreciate your advice.

>>> labelValue = 1
Loading with imageIOName: GDCM
>>> labelmapVolumeNode =slicer.mrmlScene.AddNewNodeByClass("vtkMRMLLabelMapVolumeNode")
>>> volumeNode = getNode('Displacement Field')
>>> labelmapVolumeNode = getNode('LabelMapVolume')
>>> segmentationNode = getNode('segmentation')
>>> segmentId = 'GTV'
>>> slicer.modules.segmentations.logic().ExportVisibleSegmentsToLabelmapNode(segmentationNode, labelmapVolumeNode, volumeNode)
True
>>> volumeArray = slicer.util.arrayFromVolume(volumeNode)
>>> labelArray = slicer.util.arrayFromVolume(labelmapVolumeNode)
>>> import numpy as np
>>> histogram = np.histogram(arrayFromVolume(volumeNode), bins=100)
>>> slicer.util.plot(histogram, xColumnIndex = 1)
(MRMLCorePython.vtkMRMLPlotChartNode)000001D612619E88
>>> range = [-5.0, 5.0]
>>> bins = 100
>>> histogramR, histogramBins = np.histogram(arrayFromVolume(volumeNode)[:,:,:,0], bins, range)
>>> histogramA = np.histogram(arrayFromVolume(volumeNode)[:,:,:,1], bins, range)[0]
>>> histogramS = np.histogram(arrayFromVolume(volumeNode)[:,:,:,2], bins, range)[0]
>>> slicer.util.plot(np.vstack([histogramBins[:-1], histogramR, histogramA, histogramS]).T, xColumnIndex = 0, columnNames=['N', 'R', 'A', 'S'], title='Histogram')
(MRMLCorePython.vtkMRMLPlotChartNode)000001D612619DC8
>>> magnitudes = np.linalg.norm(arrayFromVolume(volumeNode), axis=3)
>>> histogram = np.histogram(magnitudes, bins=100)
>>> slicer.util.plot(histogram, xColumnIndex = 1, columnNames=['magnitude', 'N'], title="displacement")
(MRMLCorePython.vtkMRMLPlotChartNode)000001D612619EE8
>>> magnitudes = np.linalg.norm(arrayFromVolume(volumeNode), axis=3)
>>> magnitudesInSegment = magnitudes[arrayFromSegmentBinaryLabelmap(segmentationNode, segmentId) != 0]
Traceback (most recent call last):
  File "<console>", line 1, in <module>
IndexError: boolean index did not match indexed array along dimension 0; dimension is 128 but corresponding boolean dimension is 123
>>> histogram = np.histogram(magnitudesInSegment, bins=100)
Traceback (most recent call last):
  File "<console>", line 1, in <module>
NameError: name 'magnitudesInSegment' is not defined
>>> slicer.util.plot(histogram, xColumnIndex = 1, columnNames=['magnitude', 'N'], title="displacement")
(MRMLCorePython.vtkMRMLPlotChartNode)000001D612619E28

This indicates that the size of the scalar volume and the segmentation labelmap is not the same. Using the latest Slicer Preview Release probably fixes this already, but if not then specify the displacement field as reference volume in arrayFromSegmentBinaryLabelmap.

Thank you very much.
After installing the preview version, the code also worked correctly and I was able to get the histogram successfully.

>>> labelValue = 1
Loading with imageIOName: GDCM
>>> labelmapVolumeNode =slicer.mrmlScene.AddNewNodeByClass("vtkMRMLLabelMapVolumeNode")
>>> volumeNode = getNode('Displacement Field')
>>> labelmapVolumeNode = getNode('LabelMapVolume')
>>> segmentationNode = getNode('segmentation')
>>> segmentId = 'patient'
>>> slicer.modules.segmentations.logic().ExportVisibleSegmentsToLabelmapNode(segmentationNode, labelmapVolumeNode, volumeNode)
True
>>> volumeArray = slicer.util.arrayFromVolume(volumeNode)
>>> labelArray = slicer.util.arrayFromVolume(labelmapVolumeNode)
>>> import numpy as np
>>> magnitudes = np.linalg.norm(arrayFromVolume(volumeNode), axis=3)
>>> magnitudesInSegment = magnitudes[arrayFromSegmentBinaryLabelmap(segmentationNode, segmentId) != 0]
>>> histogram = np.histogram(magnitudesInSegment, bins=100)
>>> slicer.util.plot(histogram, xColumnIndex = 1, columnNames=['magnitude', 'N'], title="displacement")
(MRMLCore.vtkMRMLPlotChartNode)000001BA99908E88
>>> segmentId = 'GTV'
>>> slicer.modules.segmentations.logic().ExportVisibleSegmentsToLabelmapNode(segmentationNode, labelmapVolumeNode, volumeNode)
True
>>> volumeArray = slicer.util.arrayFromVolume(volumeNode)
>>> labelArray = slicer.util.arrayFromVolume(labelmapVolumeNode)
>>> import numpy as np
>>> magnitudes = np.linalg.norm(arrayFromVolume(volumeNode), axis=3)
>>> magnitudesInSegment = magnitudes[arrayFromSegmentBinaryLabelmap(segmentationNode, segmentId) != 0]
>>> histogram = np.histogram(magnitudesInSegment, bins=100)
>>> slicer.util.plot(histogram, xColumnIndex = 1, columnNames=['magnitude', 'N'], title="displacement")
(MRMLCore.vtkMRMLPlotChartNode)000001BA99908C48

image

1 Like