Standard deviation of pixel values in a masked image

Question

I have a DICOM image with a mask on. It looks like a black background with a white circle in the middle (area not covered and zeroed with the mask).

The code for which is:

import numpy as np
import dicom
import pylab

ds = dicom.read_file("C:\Users\uccadmin\Desktop\James_Phantom_CT_Dec_16th\James Phantom   CT Dec 16th\Images\SEQ4Recon_3_34\IM-0268-0001.dcm")

lx, ly = ds.pixel_array.shape
X, Y = np.ogrid[0:lx, 0:ly]
mask = (X - lx/2)**2 + (Y - ly/2)**2 > lx*ly/8  # defining mask
ds.pixel_array[mask] = 0
print np.std(ds.pixel_array) # trying to get standard deviation

pylab.imshow(ds.pixel_array, cmap=pylab.cm.bone) # shows image with mask

I want to get the standard deviation of the pixel values INSIDE the white circle ONLY i.e. exclude the black space outside the circle (the mask).

I do not think the value I am getting with the above code is correct, as it is ~500, and the white circle is almost homogenous.

Any ideas how to make sure that I get the standard deviation of the pixel values within the white circle ONLY in a Pythonic way?

jwg jwg · Accepted Answer · 2014-01-24T14:09:21

I think the reason you are getting a big number is because your standard deviation is including all the zero values.

Is it enough for you to simply ignore all zero values? (This will be okay, providing that no or very few pixels in the circle have value 0.) If so

np.std([x for x in ds.pixel_array if x > 0])

should do the trick. If this isn't good enough, then you can reverse the condition in your mask to be

mask = (X - lx/2)**2 + (Y - ly/2)**2 < lx*ly/8  # defining mask, < instead of >

and do

mp.std(ds.pixel_array[mask])

Standard deviation of pixel values in a masked image

1 Answers