I can't interpret the results of a Hounsfield units conversion

[joancj97]

Context
I'm reproducing the following code from https://www.raddq.com/dicom-processing- ... in-python/. My data is a set of axial CT scans (DICOM files) from an ankle and I'm aiming to segmentate the Achilles tendon.

Code: Select all

def load_scan(path):
    slices = [pydicom.read_file(path + '/' + s) for s in os.listdir(path)]
    slices.sort(key = lambda x: int(x.InstanceNumber))
    try:
        slice_thickness = np.abs(slices[0].ImagePositionPatient[2] - slices[1].ImagePositionPatient[2])
    except:
        slice_thickness = np.abs(slices[0].SliceLocation - slices[1].SliceLocation)
        
    for s in slices:
        s.SliceThickness = slice_thickness
       
    return slices

def get_pixels_hu(scans):
    image = np.stack([s.pixel_array for s in scans])
    # Convert to int16 (from sometimes int16), 
    # should be possible as values should always be low enough (<32k)
    image = image.astype(np.int16)

    # Set outside-of-scan pixels to 1
    # The intercept is usually -1024, so air is approximately 0
    image[image == -2000] = 0
    
    # Convert to Hounsfield units (HU)
    intercept = scans[0].RescaleIntercept if 'RescaleIntercept' in scans[0] else -1024
    slope = scans[0].RescaleSlope if 'RescaleSlope' in scans[0] else 1
    
    image = slope * image.astype(np.float64)
    image = image.astype(np.int16)
    image += np.int16(intercept)
    
    return np.array(image, dtype=np.int16)

id=1
patient = load_scan(data_path)
patient_pixels = get_pixels_hu(patient)

The code uses the values Rescale Intercept (from the metadata's tag (0028,1052)) and Rescale Slope (from the metadata's tag (0028,1053)), which are absent in my dataset (it is written that the values "are stored in the DICOM header at the time of image acquisition"). I "solved" this issue by adding the if lines under ` # Convert to Hounsfield units (HU)` , but when I generate a histogram to interpret the Houndfield units values I get the following:

Code: Select all

plt.hist(patient_pixels.flatten(), bins=50, color='c')
plt.xlabel("Hounsfield Units (HU)")
plt.ylabel("Frequency")
plt.show()


The problem is that there's no tissue I found with Hounsfield units values above +1000/+3000. Thus my questions are:

Questions

1. Are the values for Rescale Intercept (-1024) and Rescale Slope (1) wrong? If so, why and which are the correct ones? If not, where I am making the mistake?
2. Does anyone know the Houndfield units for the Achilles tendon?

Other metadata values from my dataset

According to https://dicom.innolitics.com/ciods/ct-i ... e/00281052 the following information may be useful:

• Image Type (0008,0008) Value 1 is ORIGINAL and Value 3 is GDC (not LOCALIZER).
• Multi-energy CT Acquisition (0018,9361) is absent.
• If I specify Rescale Type (0028,1054) by adding

  Code: Select all

  #s.add_new([0x0028, 0x1054], 'LO', 'HU')

  in load_scan I get the same results.