In [1]:
import dicom
import pylab
import os
import numpy as np
from scipy.optimize import curve_fit
import scipy.misc
import matplotlib.pyplot as pyplot

# show dicom
def dcmshow(img):
    pylab.imshow(img.pixel_array, cmap=pylab.cm.bone)
    pylab.show()

# fitting function
def fitFunc (x, M, T1):
    return np.fabs(M * (1 - 2 * np.exp(-x/T1)))

In [2]:
rootdir = 'IR_images'

# read all files
files = []  # create an empty list
for dirName, subdirList, fileList in os.walk(rootdir):
    for filename in fileList:
        if ".dcm" in filename.lower():  # check whether the file's DICOM
            files.append(dicom.read_file(os.path.join(dirName,filename)))

In [3]:
Mimg = np.empty([files[0].Rows, files[0].Columns])
T1img = np.empty([files[0].Rows, files[0].Columns])

for r in range(0, files[0].Rows):
    for c in range(0, files[0].Columns):
        xdata = np.empty(len(files))
        ydata = np.empty(len(files))
        i = 0
        for file in files:
            xdata[i] = file.InversionTime
            ydata[i] = file.pixel_array[r][c]
            i += 1

        if(ydata[0] > 500):
            MT1 = curve_fit(fitFunc, xdata, ydata, [ydata[0], 600])[0]
            Mimg[r][c] = MT1[0]
            T1img[r][c] = MT1[1]

In [4]:
fig = pyplot.figure(dpi=256)
pyplot.set_cmap(pyplot.gray())

fig.add_subplot(121)
pyplot.imshow(Mimg)
pyplot.xlabel("Magnetization")

fig.add_subplot(122)
pyplot.xlabel("T1")
pyplot.imshow(T1img)

pyplot.show()