Convert gray pixel_array from DICOM to RGB image

I'm reading DICOM gray image file as

gray = dicom.dcmread(file).pixel_array

There I've got (x,y) shape but I need RGB (x,y,3) shape

I'm trying to convert using CV

img = cv2.cvtColor(gray, cv2.COLOR_GRAY2RGB)

And for testing I'm writing it to file cv2.imwrite('dcm.png', img)

I've got extremely dark image on output which is wrong, what is correct way to convert pydicom image to RGB?

Solution

To answer your question, you need to provide a bit more info, and be a bit clearer.

First what are you trying to do? Are you trying to only get an (x,y,3) array in memory? or are you trying to convert the dicom file to a .png file? ...they are very different things.

Secondly, what modality is your dicom image?

It's likely (unless its ultrasound or perhaps nuc med) a 16 bit greyscale image, meaning the data is 16 bit, meaning your gray array above is 16 bit data.

So the first thing to understand is window levelling and how to display a 16-bit image in 8 bits. have a look here: http://www.upstate.edu/radiology/education/rsna/intro/display.php.

If it's a 16-bit image, if you want to view it as a greyscale image in rgb format, then you need to know what window level you're using or need, and adjust appropriately before saving.

Thirdly, like lenik mention above, you need to apply the dicom slope/intercept values to your pixel data prior to using.

If your problem is just making a new array with extra dimension for rgb (so sizes (r,c) to (r,c,3)), then it's easy

# orig is your read in dcmread 2D array:
r, c = orig.shape

new = np.empty((w, h, 3), dtype=orig.dtype)

new[:,:,2] = new[:,:,1] = new[:,:,0] = orig
# or with broadcasting
new[:,:,:] = orig[:,:, np.newaxis]

That will give you the 3rd dimension. BUT the values will still all be 16-bit, not 8 bit as needed if you want it to be RGB. (Assuming your image you read with dcmread is CT, MR or equivalent 16-bit dicom - then the dtype is likely uint16).

If you want it to be RGB, then you need to convert the values to 8-bit from 16-bit. For that you'll need to decide on a window/level and apply it to select the 8-bit values from the full 16-bit data range.

Likely your problem above - I've got extremely dark image on output which is wrong - is actually correct, but it's dark because the window/level cv is using by default makes it 'look' dark, or it's correct but you didn't apply the slope/intercept.

If what you want to do is convert the dicom to png (or jpg), then you should probably use PIL or matplotlib rather than cv. Both of those offer easy ways to save a 16 bit 2D array (which is what you 'gray' is in your code above), both which allow you to specify window and level when saving to png or jpg. CV is complete overkill (meaning much bigger/slower to load, and much higher learning curve).

Some psueudo code using matplotlib. The vmin/vmax values you need to adjust - the ones here would be approximately ok for a CT image.


    import matplotlib.pyplot as plt

    df = dcmread(file)
    slope = float(df.RescaleSlope)
    intercept = float(df.RescaleIntercept)
    df_data = intercept + df.pixel_array * slope

    # tell matplotlib to 'plot' the image, with 'gray' colormap and set the
    # min/max values (ie 'black' and 'white') to correspond to 
    # values of -100 and 300 in your array
    plt.imshow(df_data, cmap='gray', vmin=-100, vmax=300)

    # save as a png file
    plt.savefig('png-copy.png')

that will save a png version, but with axes drawn as well. To save as just an image, without axes and no whitespace, use this:

    inches = (3,3)
    dpi = 150
    fig, ax = plt.subplots(figsize=inches, dpi=dpi)
    fig.subplots_adjust(left=0, right=1, top=1, bottom=0, wspace=0, hspace=0)
    ax.imshow(df_data, cmap='gray', vmin=-100, vmax=300)

    fig.save('copy-without-whitespace.png')