pythondicompydicom

Trying to input a RGB fMRI DICOM image, modify it, and save it as a grayscale dicom image in Python using pydicom

I am trying to take in an RGB fMRI scan as input and output the same scan but in grayscale with the color parts "burned" white essentially.

Whenever I try and modify any of the Data Elements, such as Photometric Interpretation and Samples Per Pixel, and use save_as to write the new DICOM file, I am unable to open that DICOM scan with the DICOM viewer giving me the error that it isn't a DICOM image.

My code is below.

import pydicom
from pydicom import dcmread
import numpy as np

#function to turn RGB array to grayscale array
#uses dot product of matrices
def rgb2gray(rgb):
   fil = [0.299, 0.587, 0.144]
   return np.dot(rgb, fil)

ds = pydicom.dcmread("dicom file")

arr = ds.pixel_array
gray_arr = rgb2gray(arr)
#gray_arr = ds.pixel_array[:,:,0]

#Have to change meta tag information when working with dicom images
ds.PhotometricInterpretation = "MONOCRHOME2"
ds.SamplesPerPixel = 1
ds.BitsAllocated = 16
ds.BitsStored = 16
ds.HighBit = 15
del ds.PlanarConfiguration
ds.is_little_endian = True
ds.fix_meta_info()


ds.PixelData = gray_arr.tobytes()
ds.save_as('fMRI.dcm', write_like_original=False)
Solution

  • The main problem is that your array has the wrong type - it is float instead of byte, so what is saved to the pixel data is a byte representation of float values (which are 4 bytes each). Also, you are setting BitsAllocated to 16, meaning that you expect 2 bytes per pixel, but in your computation you only have a byte range, e.g. you need only 8 bits per pixel.

    Lastly, you had a typo in PhotometricInterpretation ("MONOCRHOME2" instead of "MONOCHROME2"). Here is a possible solution, that converts the float array to a byte array:

    from pydicom import dcmread
from pydicom.uid import generate_uid
import numpy as np

def rgb2gray(rgb):
   fil = [0.299, 0.587, 0.144]
   return np.dot(rgb, fil)

ds = pydicom.dcmread("dicom file")

arr = ds.pixel_array
gray_arr = rgb2gray(arr).round().astype(np.uint8)

ds.PhotometricInterpretation = "MONOCHROME2"
ds.SamplesPerPixel = 1
ds.BitsAllocated = 8
ds.BitsStored = 8
ds.HighBit = 7
ds.PixelRepresenation = 0
ds.SOPInstanceUID = uid.generate_uid()
del ds.PlanarConfiguration
ds.is_little_endian = True
ds.fix_meta_info()


ds.PixelData = gray_arr.tobytes()
ds.save_as('fMRI.dcm', write_like_original=False)

    A few notes: