Is it possible to convert to Hounsfield Units in Sitk?

Beginner Questions
, ,
1

Hi,

Is there a function in SITK to convert the pixel values in dicom CT images to Hounsfield Units (HU) ? This is what I currently do:


def transform_to_hu(slices): 
    images = np.stack([file.pixel_array for file in slices]) #axis=-1 for depth in the last channel -> this yielded weired image
    images = images.astype(np.int16)

    images[images >= 1000] = 0

    for n in range(len(slices)):
        
        intercept = slices[n].RescaleIntercept
        slope = slices[n].RescaleSlope
        
        if slope != 1:
            images[n] = slope * images[n].astype(np.float64)
            images[n] = images[n].astype(np.int16)
            
        images[n] += np.int16(intercept)
    
    return np.array(images, dtype=np.int16), slices

I am specifically interested in clipping the values as follows: images[images >= 1000] = 0. When I do this on sitk to numpy converted images, it doesn’t work (plotted images are different).

2

Hello @DushiFdo,

When reading a DICOM image/series of images everything is taken care of on the ITK/SimpleITK side. The result of reading a CT will be an image with intensities in HU.

You seem to be mixing and matching SimpleITK and pydicom concepts, this is not advisable. Highly recommend that you select one tool for reading DICOM and understand its conventions.

3

Hi @zivy,

Thank you for the reply. Could you pleas tell me how to filter values in Sitk as ’ images[images >= 1000] = 0 ?

4

image*(image<1000)

1 Like
6

try
images[images > HU_max] = HU_max
images[images < HU_min] = HU_min

1 Like
7

This is a new feature in SimpleITK added for version 2.3.0 and available in the release candidates and nightly builds.