3

I am trying to make some segmentation tasks on CT images which are in DICOM format, using python.

The main problem is that when I am trying to read the DICOM file (https://www.dicomlibrary.com/meddream/md5/index.html?study=1.2.826.0.1.3680043.8.1055.1.20111102150758591.92402465.76095170) using pydicom library and use the pixel_array attribute to display the image, I`m getting an image with a much worse quality than the original (please see the code and screenshots from below).

dicom_file = pydicom.read_file('Anonymized20191102.dcm')
plt.imshow(dicom_file.pixel_array, cmap='gray')
plt.show()

The original image

The image displayed with matplotlib

My question, in order to better understand the problem is:

  1. What does pixel_array attribute represents? I have seen that it takes values in range [0, 2250].

Observations:

I also tried to get the hounsfield unit(HU) values by multiplying with slope and adding the bias, as following:

def convert_to_hu(dicom_file):
    bias = dicom_file.RescaleIntercept
    slope = dicom_file.RescaleSlope
    pixel_values = dicom_file.pixel_array
    new_pixel_values = (pixel_values * slope) + bias
    return new_pixel_values

hu_pixels = convert_to_hu(dicom_file)
plt.imshow(hu_pixels, cmap='gray')
plt.show()

The result in terms of visualization is almost the same with the one with the original values from pixel_array. Also the values are in range of [-1000, 1250] (which I think is good given the fact that HU is usually between -1000 and 1000).

robertg
  • 73
  • 1
  • 8

2 Answers2

4

The second image looks like it has plotted the pixels using a grayscale range that fits all of the provided HU values rather than just displaying the range of HU values that corresponds to the anatomy being displayed.

While I am not sure how pydicom handles this, you need to "window and level" the image to get the result from the first image. Applying a window/level, generally via a linear VOI LUT for CT images, will restrict the range of pixels displayed and scale them to the displayable range, usually 8-bit unless using grayscale monitors.

The CT image most likely has tags for Window Center (0028,1050) to set the level and Window Width (0028,1051) to set the window. The center and width values in the image were provided by the modality that generated the image and correspond to a good range to display the abdomen region in the image. The link below provides a few common window/level settings for CT images based on the body part being viewed and also provides additional information on how the settings affect the displayed image.

https://radiopaedia.org/articles/windowing-ct

While applying a window/level is useful for displaying the image, you are losing bit depth by converting from 10 or 12-bits stored per pixel (common for a CT image) to 8-bits displayed per pixel. You will need to determine if this additional bit depth is needed for your segmentation task.

2

Like Colby said, the reason it looks wrong is almost certainly the window level.

For a CT image like that, you probably want a window level around w=400, level=50, which is a min=-150, max=250. See the article Colby linked to.

You can either choose the window-level (ie min/max) values to suit, or as Colby mentioned, pull them out of the likely saved WL values in the dicom file, and calculate the min/max from them (min = level - window/2, max = level + window/2).

In your code, change to this:

level = dicom_file.WindowCenter
window = dicom_file.WindowWidth
# ...or set window/level manually to values you want
vmin = level - window/2
vmax = level + window/2
plt.imshow(hu_pixels, cmap='gray', vmin=vmin, vmax=vmax)
plt.show()

Also, you said at the start: 'I am trying to make some segmentation tasks'.

If that's the case, you need to be using the 16-bit data as-is and not an 8-bit version of, as

  1. it's got much higher data range, so is much better for segmentation, and
  2. for CT, the actual values have clinical significance, so are critical for segmentation.

...so if your overall goal is segmentation, not simply how to display the image, then displaying with the right window/level isn't really too big an issue. For segmentation you need to realise the data is full 16-bit, and use it in the full 16-bit data range.

Just figure out how WL works, and how to display it as you like/need/want with imshow() and you'll be good to go.

Amit Joshi
  • 15,448
  • 21
  • 77
  • 141
Richard
  • 3,024
  • 2
  • 17
  • 40