# Number of voxels dependent on HU value

**URL:** <https://discourse.slicer.org/t/number-of-voxels-dependent-on-hu-value/15317>\
**Category:** Support\
**Tags:** segmentation, dicom, python\
**Created:** [January 2, 2021, 8:23pm UTC](https://discourse.slicer.org/t/number-of-voxels-dependent-on-hu-value/15317 "2021-01-02T20:23:58Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![seymakandemir](https://avatars.discourse-cdn.com/v4/letter/s/a5b964/32.png) [@seymakandemir](https://discourse.slicer.org/u/seymakandemir)\
**Post date:** [January 2, 2021, 8:23pm UTC](https://discourse.slicer.org/t/number-of-voxels-dependent-on-hu-value/15317/1 "2021-01-02T20:23:58Z")

</div>

Hi all,

I created a segment of volume and got the histogram of segment. What I need to do next is to have statistical data about segment such as; How many voxels above/under 2000 HU value? And their percentage?  
I could again create 2 segments with 2000 being the boundary, but I just need the number. Is there a way that I can assing the number of voxels under 2000 HU to a variable in possibly Python Interactor.

Best regards.

Code im using:

```python
##just histogram
import numpy as np
import SampleData 

masterVolumeNode=getNode('inputVolumeNode')
segmentationNode=getNode('Segmentation')

# Create segment editor to get access to effects
segmentEditorWidget = slicer.qMRMLSegmentEditorWidget() #segment editor modulüne erişmek için
# To show segment editor widget (useful for debugging): segmentEditorWidget.show()
segmentEditorWidget.setMRMLScene(slicer.mrmlScene)
segmentEditorNode = slicer.vtkMRMLSegmentEditorNode()
slicer.mrmlScene.AddNode(segmentEditorNode)
segmentEditorWidget.setMRMLSegmentEditorNode(segmentEditorNode)
segmentEditorWidget.setSegmentationNode(segmentationNode) #segmentation
segmentEditorWidget.setMasterVolumeNode(masterVolumeNode) #segmentvolume

# Set up masking parameters
segmentEditorWidget.setActiveEffectByName("Mask volume") #mask volume kısmı histogram için segmentasyonu masklemelisin
effect = segmentEditorWidget.activeEffect()
# set fill value to be outside the valid intensity range
intensityRange = masterVolumeNode.GetImageData().GetScalarRange()
effect.setParameter("FillValue", str(intensityRange[0]-1))
# Blank out voxels that are outside the segment
effect.setParameter("Operation", "FILL_OUTSIDE")
# Create a volume that will store temporary masked volumes
maskedVolume = slicer.mrmlScene.AddNewNodeByClass("vtkMRMLScalarVolumeNode", "Temporary masked volume")
effect.self().outputVolumeSelector.setCurrentNode(maskedVolume)

# Create chart
plotChartNode = slicer.mrmlScene.AddNewNodeByClass("vtkMRMLPlotChartNode", "Histogram")

# Create histogram plot data series for each masked volume
for segmentIndex in range(segmentationNode.GetSegmentation().GetNumberOfSegments()):#for segmentIndex in segmentsayisi
  # Set active segment
  segmentID = segmentationNode.GetSegmentation().GetNthSegmentID(segmentIndex)
  segmentEditorWidget.setCurrentSegmentID(segmentID)
  # Apply mask
  effect.self().onApply()
  # Compute histogram values
  histogram = np.histogram(arrayFromVolume(maskedVolume), bins=400, range=intensityRange)
  # Save results to a new table node
  segment = segmentationNode.GetSegmentation().GetNthSegment(segmentIndex)
  tableNode=slicer.mrmlScene.AddNewNodeByClass("vtkMRMLTableNode", segment.GetName() + " histogram table")
  updateTableFromArray(tableNode, histogram)
  tableNode.GetTable().GetColumn(0).SetName("Count")
  tableNode.GetTable().GetColumn(1).SetName("Intensity")
  # Create new plot data series node
  plotSeriesNode = slicer.mrmlScene.AddNewNodeByClass("vtkMRMLPlotSeriesNode", segment.GetName() + " histogram")
  plotSeriesNode.SetAndObserveTableNodeID(tableNode.GetID())
  plotSeriesNode.SetXColumnName("Intensity")
  plotSeriesNode.SetYColumnName("Count")
  plotSeriesNode.SetPlotType(slicer.vtkMRMLPlotSeriesNode.PlotTypeScatter)
  plotSeriesNode.SetMarkerStyle(slicer.vtkMRMLPlotSeriesNode.MarkerStyleNone)
  plotSeriesNode.SetUniqueColor()
  # Add plot to chart
  plotChartNode.AddAndObservePlotSeriesNodeID(plotSeriesNode.GetID())

# Show chart in layout
slicer.modules.plots.logic().ShowChartInLayout(plotChartNode)

# Delete temporary node
slicer.mrmlScene.RemoveNode(maskedVolume)
slicer.mrmlScene.RemoveNode(segmentEditorNode)

```

 ![image](https://us1.discourse-cdn.com/flex002/uploads/slicer/original/3X/8/c/8c83f8baf75f9ac6006497b82ae15459a19518ec.png)

---

<div class="post-metadata">

**Author:** ![pll\_llq](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/pll_llq/32/9408_2.png) [@pll\_llq](https://discourse.slicer.org/u/pll_llq)\
**Post date:** [January 4, 2021, 5:49pm UTC](https://discourse.slicer.org/t/number-of-voxels-dependent-on-hu-value/15317/2 "2021-01-04T17:49:29Z")

</div>

You can convert a volume masked by a segmentation to a numpy array with `slicer.util.arrayFromVolume` and compute the number from the array.

---

<div class="post-metadata">

**Author:** ![lassoan](https://sea2.discourse-cdn.com/flex002/user_avatar/discourse.slicer.org/lassoan/32/13_2.png) [@lassoan](https://discourse.slicer.org/u/lassoan)\
**Post date:** [January 4, 2021, 10:20pm UTC](https://discourse.slicer.org/t/number-of-voxels-dependent-on-hu-value/15317/3 "2021-01-04T22:20:44Z")

</div>

You have already solved the problem. You can get the histogram with any resolution and then sum values of the histogram to get total number of voxels in a certain intensity range.

Alternatively, you can find a complete implementation of this feature as a module with a nice GUI in [LungCTAnalyzer extension](https://discourse.slicer.org/t/new-lungctanalyzer-extension-for-lung-ct-segmentation-and-analysis-for-covid-19-assessment/15006): it computes volume of each region that is specified by voxel intensity range. It shows nice real-time preview, generates, PDF report, etc. To get started, you can use the module as is, just ignoring that takes “left lung” and “right lung” segments as inputs. You can put any other structure in those segments and you can specify any intensity range. If it is close enough to what you want to achieve, you can easily customize it to be more convenient and intuitive for your data sets - it is just a Python scripted module that you can modify with a text editor.

---

<div class="post-metadata">

**Author:** ![seymakandemir](https://avatars.discourse-cdn.com/v4/letter/s/a5b964/32.png) [@seymakandemir](https://discourse.slicer.org/u/seymakandemir)\
**Post date:** [January 5, 2021, 10:30am UTC](https://discourse.slicer.org/t/number-of-voxels-dependent-on-hu-value/15317/5 "2021-01-05T10:30:14Z")

</div>

> [@seymakandemir](#):
>
> I looked at the LungCT analyzer you mentioned. Using this module, I can get the maximum HU value of 1000. I could not use it because I wanted to work in a wider range.  
> I changed the code to get the voxel number above and below 2000HU. I still couldn’t make it look like this percentage is below or above 2000HU.  
> My modified code:

import numpy as np  
import SampleData  
masterVolumeNode=getNode(‘inputVolumeNode\_masked’)  
segmentationNode=getNode(‘Segmentation’)

# Create segment editor to get access to effects

segmentEditorWidget = slicer.qMRMLSegmentEditorWidget() #segment editor modulüne erişmek için

# To show segment editor widget (useful for debugging): segmentEditorWidget.show()

segmentEditorWidget.setMRMLScene(slicer.mrmlScene)  
segmentEditorNode = slicer.vtkMRMLSegmentEditorNode()  
slicer.mrmlScene.AddNode(segmentEditorNode)  
segmentEditorWidget.setMRMLSegmentEditorNode(segmentEditorNode)  
segmentEditorWidget.setSegmentationNode(segmentationNode) #segmentation  
segmentEditorWidget.setMasterVolumeNode(masterVolumeNode) #segmentvolume

# Set up masking parameters

segmentEditorWidget.setActiveEffectByName(“Mask volume”) #mask volume kısmı histogram için segmentasyonu masklemelisin  
effect = segmentEditorWidget.activeEffect()

# set fill value to be outside the valid intensity range

intensityRange = masterVolumeNode.GetImageData().GetScalarRange()  
effect.setParameter(“FillValue”, str(intensityRange[0]-1))

# Blank out voxels that are outside the segment

effect.setParameter(“Operation”, “FILL\_OUTSIDE”)

# Create a volume that will store temporary masked volumes

maskedVolume = slicer.mrmlScene.AddNewNodeByClass(“vtkMRMLScalarVolumeNode”, “Temporary masked volume”)  
effect.self().outputVolumeSelector.setCurrentNode(maskedVolume)

# Create chart

#plotChartNode = slicer.mrmlScene.AddNewNodeByClass(“vtkMRMLPlotChartNode”, “Histogram”)

# Create histogram plot data series for each masked volume

for segmentIndex in range(segmentationNode.GetSegmentation().GetNumberOfSegments()):#for segmentIndex in segmentsayisi

# Set active segment

segmentID = segmentationNode.GetSegmentation().GetNthSegmentID(segmentIndex)  
segmentEditorWidget.setCurrentSegmentID(segmentID)

# Apply mask

effect.self().onApply()

# Compute histogram values

histogram\_upper = np.histogram(arrayFromVolume(maskedVolume), bins=1, range=(2000,np.amax(intensityRange)))

# Save results to a new table node

segment = segmentationNode.GetSegmentation().GetNthSegment(segmentIndex)  
tableNode=slicer.mrmlScene.AddNewNodeByClass(“vtkMRMLTableNode”, segment.GetName() + " histogram\_upper table")  
updateTableFromArray(tableNode, histogram\_upper)  
tableNode.GetTable().GetColumn(0).SetName(“Count”)  
tableNode.GetTable().GetColumn(1).SetName(“Intensity”)

# Create new plot data series node

#histogram\_lower  
histogram\_lower= np.histogram(arrayFromVolume(maskedVolume), bins=1, range=(np.amin(intensityRange),1999))  
tableNode=slicer.mrmlScene.AddNewNodeByClass(“vtkMRMLTableNode”, segment.GetName() + " histogram\_lower table")  
updateTableFromArray(tableNode, histogram\_lower)  
tableNode.GetTable().GetColumn(0).SetName(“Count”)  
tableNode.GetTable().GetColumn(1).SetName(“Intensity”)

# Create new plot data series node

# Delete temporary node

slicer.mrmlScene.RemoveNode(maskedVolume)  
slicer.mrmlScene.RemoveNode(segmentEditorNode)

Fig.1: Number of voxels under 2000HU

 ![image](https://us1.discourse-cdn.com/flex002/uploads/slicer/original/3X/6/0/60a3821227e256a9faba0b3610b1b6d9827443b0.png)  
Fig.2: Number of voxels above 2000HU  
 ![image](https://us1.discourse-cdn.com/flex002/uploads/slicer/original/3X/8/3/83cf31903f55cc8ccd035bd4996de0e9ac41d653.png)

Thank you

---

<div class="post-metadata">

**Author:** ![seymakandemir](https://avatars.discourse-cdn.com/v4/letter/s/a5b964/32.png) [@seymakandemir](https://discourse.slicer.org/u/seymakandemir)\
**Post date:** [January 5, 2021, 10:33am UTC](https://discourse.slicer.org/t/number-of-voxels-dependent-on-hu-value/15317/6 "2021-01-05T10:33:00Z")

</div>

Since I have just started using this program, I am not very familiar yet. I don’t know how to use this.  
Thank you