HW5 Prob 7.4 Python plotting code

hw5_p7_4_plot.py

# Plotting Code
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.colors
import scipy.io
plt.style.use('bmh')

k = 50 # image size is k x k
m = k * k # total number of pixels

# define regions
tumor   = np.full([k, k], False, dtype=np.bool)
crit    = np.full([k, k], False, dtype=np.bool)

tumor[21:30, 20:32] = True    # tumor
crit[33:40, 16:25] = True     # critical structure 1
crit[23:32, 38:45] = True     # critical structure 2
other = np.logical_not(np.logical_or(tumor, crit)) # other region

tumor_min = 0.9; tumor_max = 1.5
other_min = 0.0; other_max = 0.5 
crit_min  = 0.0; crit_max  = 0.1 

Dmin = np.zeros([k, k])
Dmin[other] = other_min    # regular cells
Dmin[tumor] = tumor_min    # tumor
Dmin[crit]  = crit_min     # critical regions

Dmax = np.zeros([k, k]);
Dmax[other] = other_max    # regular cells
Dmax[tumor] = tumor_max    # tumor
Dmax[crit]  = crit_max     # critical regions

Ntheta = 20   # number of angles (equally spaced from 0 to pi)
Nd = 20       # number of parallel lines for each angle
n = Nd*Ntheta # total number of measurements

# Load dosage data from file
dosage_data = scipy.io.loadmat('dosage_opt_data.mat')
A = dosage_data['A']

# dmin and dmax contain the vectorized version of the images Dmin and Dmax
dmin = Dmin.flatten()
dmax = Dmax.flatten()

# ----------------------------------------------------------
# Plot area geometry
# ----------------------------------------------------------
# Sets the colormap so that:
#   tumor regions are in red
#   critical regions are in blue
#   other regions are in green
C = matplotlib.colors.LinearSegmentedColormap.from_list(name='C', colors=[
    (0.0, 0.0, 0.5),
    (0.0, 0.5, 0.0),
    (0.0, 0.5, 0.0),
    (0.0, 0.5, 0.0),
    (0.5, 0.0, 0.0)])
f, ax = plt.subplots(1, 1)
ax.imshow(Dmin, cmap=C)
ax.set_title('Area Geometry')
f.savefig('Dmax.png')

#----------------------------------------------------------
# Plot dosage level for uniform beam levels
# ----------------------------------------------------------
b = 0.05 * np.ones([n,1]) # uniform level of 0.05 for all beams
d = np.matmul(A, b)
D = d.reshape([k, k])

f, ax = plt.subplots(1, 1)
ax.imshow(D, cmap=C)
ax.set_title('Dosage (uniform beam pattern)')
f.savefig('unif_D.png')

# ----------------------------------------------------------
# Plot histogram of pixel dosages for the three regions
# ----------------------------------------------------------
# You may simply copy and paste the following code to plot 
# the histogram for your solution (assuming you store it in
# the variable d).
#
f, axs = plt.subplots(3,1)
f.suptitle('Dosage histograms (uniform beam pattern)') # set title of figure

axs[0].hist(D[tumor], bins=np.arange(0, 1.5, 0.05), color=[0.5, 0, 0])
axs[0].axvline(tumor_min - 0.025, linewidth=2, c='k', ls='dashed')
axs[0].axvline(tumor_max + 0.025, linewidth=2, c='k', ls='dashed')
axs[0].set_xlim([-0.1, 1.6])
axs[0].set_ylabel('tumor')

axs[1].hist(D[crit], bins=np.arange(0, 1.5, 0.05), color=[0, 0, 0.5])
axs[1].axvline(crit_min - 0.025, linewidth=2, c='k', ls='dashed')
axs[1].axvline(crit_max + 0.025, linewidth=2, c='k', ls='dashed')
axs[1].set_xlim([-0.1, 1.6])
axs[1].set_ylabel('critical')

axs[2].hist(D[other], bins=np.arange(0, 1.5, 0.05), color=[0, 0.5, 0])
axs[2].axvline(other_min - 0.025, linewidth=2, c='k', ls='dashed')
axs[2].axvline(other_max + 0.025, linewidth=2, c='k', ls='dashed')
axs[2].set_xlim([-0.1, 1.6])
axs[2].set_ylabel('other')
f.savefig('unif_hist.png')

Tested with

Python:  3.6.5 (default, Mar 30 2018, 06:41:53) 
[GCC 4.2.1 Compatible Apple LLVM 9.0.0 (clang-900.0.39.2)]
numpy:  1.14.3
matplotlib:  2.2.2
scipy:  1.1.0rc1