<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Lanex Strip Test
Some quick code for analyzing the data acquired during the Lanex strip test. In this test, the lanex strip used in the electron spectrometer design was partially covered with pieces of electrical tape before being exposed to the electron source. This was done to see how much of the signal is due to Lanex luminesence and how much is due to direct excitation in the CCD by electrons.
<img src="lanex_strips.jpg" height="628" width="368" />
Step5:  Espec functions
Step11:  Data read functions
Step12:  Load data and subtract background
Step13:  Mesh plot all unsaturated data
Step14:  Plot averaged data
Step15:  Interpolate data between regions with and without vinyl
Step16:  Plot interpolated data
Step17:  Calculate signal ratios
Step18:  Determine signal ratios from MCNP output
<font size="4"><p>MCNP was ran at the above chosen energies both with and without vinyl covering the lanex. The CCD depletion region thickness was also varied from 1 um to 20 um for each energy. The signal ratio at each of these energies can be determined by finding the ratio between the total signal without vinyl and the total signal with vinyl. This can then be compared with the experimentally determined signal ratio in order to find the appropriate depletion region thickness.</p>
<p>First, we need to read in the MCNP data for the "without vinyl" case. The energy deposited into each segment of the phosphor layer can be added to an array (code adapted from "Read MCNP Output" script)
Step19:  <font size="4"><p>Plot energy absorbed versus segment depth to make sure everything looks good
Step20:  <font size="4"><p>Next, the number of photons that reach the CCD can be calculated. From the quantum efficiency of the CCD, the contribution to the signal from visible photons can be found.</p>
Step21:  <font size="4"><p>The amount of energy deposited into the CCD depletion region from electrons and x-rays can be read for both the "with vinyl" and "without vinyl" cases
Step22:  <font size="4"><p>Finally, determine the signal ratios</p>
Step23:  <font size="4"><p>and plot the results
<ASSISTANT_TASK:>
Python Code:
#Imports
from math import *
import numpy as np
import scipy as sp
import scipy.special
import scipy.interpolate as interpolate
import matplotlib as mpl
import matplotlib.pyplot as plt
import matplotlib.cbook as cbook
import seaborn as sns
import sys
import os
#Import custom modules
from physics import *
%matplotlib notebook
Explanation: Lanex Strip Test
Some quick code for analyzing the data acquired during the Lanex strip test. In this test, the lanex strip used in the electron spectrometer design was partially covered with pieces of electrical tape before being exposed to the electron source. This was done to see how much of the signal is due to Lanex luminesence and how much is due to direct excitation in the CCD by electrons.
<img src="lanex_strips.jpg" height="628" width="368" />
End of explanation
B0 = 1710.0/10**4 #Magnetic field strength in Tesla
def KEcalc(z,y):
    Returns KE in J given z-position in m
    return me*c**2*(sqrt((q*B0/(me*c))**2*((z**2+y**2)/(2*y))**2+1)-1)
def Radius(KE):
    #Radius of electron orbit in m given KE in keV
    return me*c/(q*B0)*sqrt((KE*1000*q/(me*c**2)+1)**2-1)
def zfcalc(KE,y):
    #Returns z-position at screen in inches given KE in keV
    R = Radius(KE)
    return sqrt(R**2 - (y-R)**2)
def zfcalcGeneral(KE,yM,y):
    R = Radius(KE)
    zM = zfcalc(KE,yM)
    return zM + (y - yM)*(R - yM)/zM
def KEcalcGeneral(zf,yM,yf):
    Returns KE in J given z-position of electrons, y-position of magnet edge, and y-position of screen, all in m
    a = (yM+yf)**2
    b = -2*yM*(yf*(yM+yf)+zf**2)
    d = yM**2*(zf**2+yf**2)
    f = (me*c)/(q*B0)
    g = (-b+sqrt(b**2-4*a*d))/(2*a)
    return me*c**2*(sqrt(g**2+f**2)/f - 1)
def AngleIncidence(KE,yM):
    R = Radius(KE)
    return asin((R-yM)/R)
Explanation: Espec functions
End of explanation
def getfns(folder,ext=''):
    Get a list of full path filenames for all files in a folder and subfolders for the given extension
    fns = []
    for file in os.listdir(folder):
            if file.endswith(ext):
                fns.append(os.path.join(folder,file))
    return fns
def readcsv(filename):
    Read in a csv file and load it as an array
    return np.loadtxt(open(filename, "rb"), delimiter=",")
def readcsvs(filenames):
    Read in multiple csv files and load them in to a two-dimensional array
    template = readcsv(filenames[0])
    numfiles = len(filenames)
    Data = np.zeros([numfiles,len(template)])
    for i in range(numfiles):
        spectrum = readcsv(filenames[i])
        Data[i,:] = spectrum
    return Data
def DataClean(Data):
    Read in data and clean it by removing rows that have saturated pixel values
    maxes = np.max(Data[:,500:],1)
    includes = maxes<(2**16-1)
    rejects = maxes>(2**16-2)
    CleanData = Data[includes,:]
    return CleanData
def DataAverage(Data):
    Average input 2d array into a 1d array
    return np.mean(Data,0)
Explanation: Data read functions
End of explanation
MagnetOutFolderPath = os.curdir + '/Data/2015-08-17_Tiger_stripes_test/even_more_no_magnet'
MagnetInFolderPath = os.curdir + '/Data/2015-08-17_Tiger_stripes_test/magnet_in'
NoLaserFolderPath = os.curdir + '/Data/2015-08-17_Tiger_stripes_test/no_laser'
MagnetOutFiles = getfns(MagnetOutFolderPath,'')
MagnetInFiles = getfns(MagnetInFolderPath,'')
NoLaserFiles = getfns(NoLaserFolderPath,'')
MagnetOutData = readcsvs(MagnetOutFiles)
MagnetInData = readcsvs(MagnetInFiles)
Background = DataAverage(MagnetInData)
MagnetOutDataClean = DataClean(MagnetOutData)-Background
Explanation: Load data and subtract background
End of explanation
sns.set(font_scale=1.5)
fig1 = plt.figure(figsize=(8,6))
ax1 = fig1.add_subplot(111)
ax1.set_xlim(500,3648)
ax1.set_ylim(0,len(MagnetInData))
ax1.set_xlabel('Pixel')
ax1.set_ylabel('Instance')
# mesh1 = ax1.pcolormesh(MagnetInData, cmap='hot',vmin=0, vmax=10000)
mesh1 = ax1.pcolormesh(MagnetOutDataClean, cmap='inferno',vmin=0, vmax=65535)
Explanation: Mesh plot all unsaturated data
End of explanation
MagnetOutDataCleanAvg = DataAverage(MagnetOutDataClean)
MagnetOutDataCleanAvgHat = savitzky_golay(MagnetOutDataCleanAvg,51,3) #Smoothed Data
fig2 = plt.figure(figsize=(12,8))
ax2 = fig2.add_subplot(111)
ax2.set_xlim(0,3648)
ax2.set_ylim(0,65535)
ax2.set_xlabel('Pixel Number')
ax2.set_ylabel('Pixel Value')
#ax2.semilogy()
ax2.plot(MagnetOutDataCleanAvg,linewidth=1)
ax2.plot(MagnetOutDataCleanAvgHat,linewidth=1,color='r')
Explanation: Plot averaged data
End of explanation
#Pixel values were determined by picking a region within each peak and trough
LanexPixelArray = np.hstack((np.arange(240,260),np.arange(670,750),np.arange(1090,1300),np.arange(1600,1700),\
                            np.arange(1930,2030),np.arange(2400,2450),np.arange(2820,2920),np.arange(3280,3350)))
LanexData = np.hstack((MagnetOutDataCleanAvgHat[240:260],MagnetOutDataCleanAvgHat[670:750],\
                       MagnetOutDataCleanAvgHat[1090:1300],MagnetOutDataCleanAvgHat[1600:1700],\
                       MagnetOutDataCleanAvgHat[1930:2030],MagnetOutDataCleanAvgHat[2400:2450],\
                       MagnetOutDataCleanAvgHat[2820:2920],MagnetOutDataCleanAvgHat[3280:3350]))
VinylPixelArray = np.hstack((np.arange(460,490),np.arange(900,950),np.arange(1400,1500),np.arange(1750,1850),\
                            np.arange(2150,2250),np.arange(2540,2640)))
VinylData = np.hstack((MagnetOutDataCleanAvgHat[460:490],MagnetOutDataCleanAvgHat[900:950],\
                       MagnetOutDataCleanAvgHat[1400:1500],MagnetOutDataCleanAvgHat[1750:1850],\
                       MagnetOutDataCleanAvgHat[2150:2250],MagnetOutDataCleanAvgHat[2540:2640]))
LanexInterpFunc = interpolate.interp1d(LanexPixelArray,LanexData,kind='slinear')
LanexPixelArrayInterp = np.arange(240,3350)
LanexDataInterp = LanexInterpFunc(LanexPixelArrayInterp)
VinylInterpFunc = interpolate.interp1d(VinylPixelArray,VinylData,kind='slinear')
VinylPixelArrayInterp = np.arange(460,2640)
VinylDataInterp = VinylInterpFunc(VinylPixelArrayInterp)
Explanation: Interpolate data between regions with and without vinyl
End of explanation
# sns.set(context='poster',font_scale=1.5)
# sns.set_style("darkgrid")
# sns.set_palette(palette='deep')
# sns.set_color_codes(palette='deep')
plt.figure(figsize=(12,6))
#With Lanex
#plt.plot(LanexPixelArray,LanexData,linestyle='None',marker='.')
plt.plot(LanexPixelArrayInterp,LanexDataInterp, linewidth=2,linestyle='--',label='Without Vinyl',color='b')
#Without Lanex
#plt.plot(VinylPixelArray,VinylData,linestyle='None',marker='.')
plt.plot(VinylPixelArrayInterp,VinylDataInterp,linewidth=2,linestyle='--',label='With Vinyl',color='r')
#All smoothed data
plt.plot(MagnetOutDataCleanAvgHat,linewidth=2,label='All Data',color='g')
plt.xlim(0,2000)
plt.xlabel('Pixel Number')
plt.ticklabel_format(style='sci', axis='y', scilimits=(0,0))
plt.ylabel('Pixel Value')
plt.legend()
plt.subplots_adjust(left=0.12,bottom=0.14) #Adjust spacing to prevent clipping of x and y labels
plt.savefig('strip_test.svg')
Explanation: Plot interpolated data
End of explanation
yM = 0.5 #magnet edge position in inches
CCDpos = 3.02 #CCD y-position relative to magnet edge in mm
yM = yM*.0254 #Convert to base units
CCDpos = CCDpos*10**-3 #Convert to base units
yf = yM + CCDpos #Screen position
KEcalcGeneralVec = np.vectorize(KEcalcGeneral) #Vectorize espec function for KE
ChosenPixels = [479.37419602, 658.07909036, 805.53608568, 929.62974462, 1079.52432829, 1325.661891,\
                1454.37508809, 1637.97442378, 1801.34325342, 1984.12923212, 2202.34572023, 2412.63089924,\
                2595.3896852]
ChosenPositions = np.multiply(np.add(np.multiply(ChosenPixels,8.0/10**3),9.5),10**-3) #Convert to meters
ChosenEnergies = KEcalcGeneralVec(ChosenPositions,yM,yf) #Energies corresponding to the chosen pixels
SignalRatio = LanexInterpFunc(ChosenPixels)/VinylInterpFunc(ChosenPixels)
print(ChosenEnergies/(1000*q),'\n')
print(SignalRatio)
Explanation: Calculate signal ratios
End of explanation
mass = 2.95680*10**-2 #mass of active layer in g
Nsegments = 101 #number of segments
segmentmass = mass/Nsegments #mass of each segment in g
phosphorthickness = 81.*10**-6 #thickness of active layer in m
KEsim = np.multiply([338, 391, 440, 485],10**-3) #Chosen energies (only the first 4 were used to save time)
#Create array of angles
CosThetaSim = []
for x in KEsim:
    CosThetaSim.append(np.cos(AngleIncidence(x*1000,yM)))
    
CCDthickness = np.arange(1,21,1) #Depletion region thickness array
MCNP_Directory = '/home/drake/Documents/Physics/Research/Python/MCNP_Code'
#Read files. As far as visible photons are concerned, the depletion region thickness doesn't matter
deptharray= []
Eabsorbed_Segments = []
Eabsorbed_Segments_Error = []
for i in range(len(KEsim)):
    deptharray.append([])
    Eabsorbed_Segments.append([])
    Eabsorbed_Segments_Error.append([])
    #Directory of the input files where vinyl was not used in the simulation
    
    if i%2==0: #Energies where vinyl was experimentally used
        directory = MCNP_Directory + '/MCNP_Decks_Varying_CCD_Thickness/Old_Output/Output2/Out_{KE}MeV_{Theta}Degrees_20umCCD_inverted'\
        .format(KE=str(round(KEsim[i],3)),Theta=str(round(np.arccos(CosThetaSim[i])*360/(2*pi),1)))
        
    else: #Energies where vinyl wasn't experimentally used
        directory = MCNP_Directory + '/MCNP_Decks_Varying_CCD_Thickness/Old_Output/Output1/Out_{KE}MeV_{Theta}Degrees_20umCCD'\
        .format(KE=str(round(KEsim[i],3)),Theta=str(round(np.arccos(CosThetaSim[i])*360/(2*pi),1)))
    for segmentnumber in range(Nsegments-1):
        segment = segmentnumber+100 #segment number
        deptharray[i].append(segmentnumber*phosphorthickness/Nsegments) #depth for each segment in microns
        printflag = False
        with open(directory) as searchfile:
            for line in searchfile:
                left,sep,right = line.partition('       -{segmentlabel}  '.format(segmentlabel=str(segment)))
                if printflag:
                    Eabs_per_g = float(line[17:28])
                    Eabs_per_g_Error = float(line[29:35])
                    Eabsorbed_Segments[i].append(Eabs_per_g*segmentmass)
                    Eabsorbed_Segments_Error[i].append(Eabs_per_g_Error*segmentmass)
                    printflag = False
                if sep:
                    printflag = True
Explanation: Determine signal ratios from MCNP output
<font size="4"><p>MCNP was ran at the above chosen energies both with and without vinyl covering the lanex. The CCD depletion region thickness was also varied from 1 um to 20 um for each energy. The signal ratio at each of these energies can be determined by finding the ratio between the total signal without vinyl and the total signal with vinyl. This can then be compared with the experimentally determined signal ratio in order to find the appropriate depletion region thickness.</p>
<p>First, we need to read in the MCNP data for the "without vinyl" case. The energy deposited into each segment of the phosphor layer can be added to an array (code adapted from "Read MCNP Output" script):</p>
</font>
End of explanation
plt.figure(figsize=(12,6))
plt.plot(np.multiply(deptharray[1],10**6),np.multiply(Eabsorbed_Segments[0],10**3),label='338 keV')
plt.plot(np.multiply(deptharray[1],10**6),np.multiply(Eabsorbed_Segments[1],10**3),label='391 keV')
plt.plot(np.multiply(deptharray[1],10**6),np.multiply(Eabsorbed_Segments[2],10**3),label='440 keV')
plt.plot(np.multiply(deptharray[1],10**6),np.multiply(Eabsorbed_Segments[3],10**3),label='485 keV')
plt.xlabel('Depth (um)')
plt.ylabel('Energy absorbed (keV)')
plt.legend()
Explanation: <font size="4"><p>Plot energy absorbed versus segment depth to make sure everything looks good:</p></font>
End of explanation
scatteringlength = 2.84*10**-6 #photon scattering length in Gd2O2S in m
conversionefficiency = 0.16 #electron energy to light energy conversion efficiency
emissionwavelength = 545*10**-9 #lanex emission wavelength in m
emissionenergy = hbar*2*pi*c/emissionwavelength
photonNumberArray = []
for i in range(len(KEsim)):
    photonNumber = 0
    for j in range(len(deptharray[i])):
        photonNumber += conversionefficiency*Eabsorbed_Segments[i][j]*10**6*q/emissionenergy\
        *(j+0.5)/Nsegments # Nabs = Nexc*(Distance from top of lanex)/(Phosphor thickness)
        #where Distance from top of lanex = (SegmentNumber+0.5)/(Nsegments)*(Phosphor thickness)
    photonNumberArray.append(photonNumber)
QuantumEfficiency = 0.4
PhotonSignal = np.array(np.multiply(photonNumberArray,QuantumEfficiency))
print(PhotonSignal)
Explanation: <font size="4"><p>Next, the number of photons that reach the CCD can be calculated. From the quantum efficiency of the CCD, the contribution to the signal from visible photons can be found.</p>
End of explanation
#With vinyl
ExAbsorbedCCDVinyl = [] #x-rays
ExAbsorbedCCDErrorVinyl = []
EelAbsorbedCCDVinyl = [] #electrons
EelAbsorbedCCDErrorVinyl = []
CCDmass = np.multiply(CCDthickness,1.39740*10**-4) #mass of photoactive layer of CCD
for i in range(len(KEsim)):
    
    ExAbsorbedCCDVinyl.append([])
    ExAbsorbedCCDErrorVinyl.append([])
    EelAbsorbedCCDVinyl.append([])
    EelAbsorbedCCDErrorVinyl.append([])
    
    for j in range(len(CCDthickness)):
    
        #Directory of the input files where vinyl was used in the simulation
        
        if i%2==0: #Energies where vinyl was experimentally used
            directory = MCNP_Directory + '/MCNP_Decks_Varying_CCD_Thickness/Old_Output/Output1/Out_{KE}MeV_{Theta}Degrees_{CCD}umCCD'\
            .format(KE=str(round(KEsim[i],3)),Theta=str(round(np.arccos(CosThetaSim[i])*360/(2*pi),1)),\
                   CCD=str(CCDthickness[j]))
        else: #Energies where vinyl wasn't experimentally used
            directory = MCNP_Directory + '/MCNP_Decks_Varying_CCD_Thickness/Old_Output/Output2/Out_{KE}MeV_{Theta}Degrees_{CCD}umCCD_inverted'\
            .format(KE=str(round(KEsim[i],3)),Theta=str(round(np.arccos(CosThetaSim[i])*360/(2*pi),1)),\
                   CCD=str(CCDthickness[j]))
        printflag = False
        with open(directory) as searchfile:
            firstoccurence = False
            for line in searchfile:
                left,sep,right = line.partition('cell  7')
                if printflag:
                    Eabs_per_g = float(line[17:28])
                    Eabs_per_g_Error = float(line[29:35])
                    if firstoccurence:
                        ExAbsorbedCCDVinyl[i].append(Eabs_per_g*CCDmass[j])
                        ExAbsorbedCCDErrorVinyl[i].append(Eabs_per_g_Error*CCDmass[j])
                    else:
                        EelAbsorbedCCDVinyl[i].append(Eabs_per_g*CCDmass[j])
                        EelAbsorbedCCDErrorVinyl[i].append(Eabs_per_g_Error*CCDmass[j])
                    printflag = False
                if sep: # True iff 'cell  6' in line
                    printflag = True
                    firstoccurence = not firstoccurence
                    
#Calculate the overall contribution to the signal
Energy_eh_pair = 3.65 #Energy needed to generate electron-hole pair in Si in eV
ElectronSignalVinyl = []
XraySignalVinyl = []
for i in range(len(KEsim)):
    ElectronSignalVinyl.append([])
    XraySignalVinyl.append([])
    
    ElectronSignalVinyl[i].append(np.multiply(EelAbsorbedCCDVinyl[i],10**6/Energy_eh_pair))
    XraySignalVinyl[i].append(np.multiply(ExAbsorbedCCDVinyl[i],10**6/Energy_eh_pair))
ElectronSignalVinyl = np.array(ElectronSignalVinyl)
XraySignalVinyl = np.array(XraySignalVinyl)
TotalSignalVinyl = np.add(ElectronSignalVinyl,XraySignalVinyl)
#Without vinyl
ExAbsorbedCCDNoVinyl = [] #x-rays
ExAbsorbedCCDErrorNoVinyl = []
EelAbsorbedCCDNoVinyl = [] #electrons
EelAbsorbedCCDErrorNoVinyl = []
CCDmass = np.multiply(CCDthickness,1.39740*10**-4) #mass of photoactive layer of CCD
for i in range(len(KEsim)):
    
    ExAbsorbedCCDNoVinyl.append([])
    ExAbsorbedCCDErrorNoVinyl.append([])
    EelAbsorbedCCDNoVinyl.append([])
    EelAbsorbedCCDErrorNoVinyl.append([])
    
    for j in range(len(CCDthickness)):
    
        #Directory of the input files where vinyl was not used in the simulation
        
        if i%2==0: #Energies where vinyl was experimentally used
            directory = MCNP_Directory + '/MCNP_Decks_Varying_CCD_Thickness/Old_Output/Output2/Out_{KE}MeV_{Theta}Degrees_{CCD}umCCD_inverted'\
            .format(KE=str(round(KEsim[i],3)),Theta=str(round(np.arccos(CosThetaSim[i])*360/(2*pi),1)),\
                   CCD=str(CCDthickness[j]))
        else: #Energies where vinyl wasn't experimentally used
            directory = MCNP_Directory + '/MCNP_Decks_Varying_CCD_Thickness/Old_Output/Output1/Out_{KE}MeV_{Theta}Degrees_{CCD}umCCD'\
            .format(KE=str(round(KEsim[i],3)),Theta=str(round(np.arccos(CosThetaSim[i])*360/(2*pi),1)),\
                   CCD=str(CCDthickness[j]))
        printflag = False
        with open(directory) as searchfile:
            firstoccurence = False
            for line in searchfile:
                left,sep,right = line.partition('cell  6')
                if printflag:
                    Eabs_per_g = float(line[17:28])
                    Eabs_per_g_Error = float(line[29:35])
                    if firstoccurence:
                        ExAbsorbedCCDNoVinyl[i].append(Eabs_per_g*CCDmass[j])
                        ExAbsorbedCCDErrorNoVinyl[i].append(Eabs_per_g_Error*CCDmass[j])
                    else:
                        EelAbsorbedCCDNoVinyl[i].append(Eabs_per_g*CCDmass[j])
                        EelAbsorbedCCDErrorNoVinyl[i].append(Eabs_per_g_Error*CCDmass[j])
                    printflag = False
                if sep: # True iff 'cell  6' in line
                    printflag = True
                    firstoccurence = not firstoccurence
                    
#Calculate the overall contribution to the signal
Energy_eh_pair = 3.65 #Energy needed to generate electron-hole pair in Si in eV
ElectronSignalNoVinyl = []
XraySignalNoVinyl = []
for i in range(len(KEsim)):
    ElectronSignalNoVinyl.append([])
    XraySignalNoVinyl.append([])
    
    ElectronSignalNoVinyl[i].append(np.multiply(EelAbsorbedCCDNoVinyl[i],10**6/Energy_eh_pair))
    XraySignalNoVinyl[i].append(np.multiply(ExAbsorbedCCDNoVinyl[i],10**6/Energy_eh_pair))
    
TotalSignalNoVinyl = np.add(ElectronSignalNoVinyl,XraySignalNoVinyl)
for i in range(len(KEsim)):
    TotalSignalNoVinyl[i] += PhotonSignal[i]
Explanation: <font size="4"><p>The amount of energy deposited into the CCD depletion region from electrons and x-rays can be read for both the "with vinyl" and "without vinyl" cases:</p>
End of explanation
Ratio = [TotalSignalNoVinyl[0][0][0:19]/TotalSignalVinyl[0][0][0:19],\
                      TotalSignalNoVinyl[1][0][0:19]/TotalSignalVinyl[1][0][0:19],\
                      TotalSignalNoVinyl[2][0][0:19]/TotalSignalVinyl[2][0][0:19],\
                      TotalSignalNoVinyl[3][0][0:19]/TotalSignalVinyl[3][0][0:19]]
MCNPSignalRatio = []
for i in range(len(CCDthickness)-1):
    MCNPSignalRatio.append([])
    for j in range(len(KEsim)):
        MCNPSignalRatio[i].append(Ratio[j][i])
Explanation: <font size="4"><p>Finally, determine the signal ratios</p>
End of explanation
mpl.rcParams.update({'font.size': 24, 'font.family': 'serif'})
sns.set(context='poster',font_scale=1.5)
sns.set_style("darkgrid")
sns.set_palette(palette='deep')
sns.set_color_codes(palette='deep')
plt.figure(figsize=(12,6))
plt.plot(np.multiply(KEsim,1.0),MCNPSignalRatio[0],linewidth=2,label='1 '+ u'\u03bcm')
plt.plot(np.multiply(KEsim,1.0),MCNPSignalRatio[1],linewidth=2,label='2 '+ u'\u03bcm')
plt.plot(np.multiply(KEsim,1.0),MCNPSignalRatio[4],linewidth=2,color='y',label='5 '+ u'\u03bcm')
plt.plot(np.multiply(KEsim,1.0),MCNPSignalRatio[9],linewidth=2,color='c',label='10 '+ u'\u03bcm')
#plt.plot(np.multiply(KEsim,1.0),MCNPSignalRatio[14],linewidth=2,color='purple',label='MCNP 15 '+ u'\u03bcm photoactive region')
plt.plot(np.divide(ChosenEnergies,10**6*q),SignalRatio,linewidth=2,color='r',marker='o',label='Experiment')
plt.xlabel('Electron Energy (keV)')
plt.ylabel('Signal Ratio (No Vinyl / Vinyl)')
plt.xlim(0.33,0.50)
plt.legend(title='Photoactive Region Thickness')
plt.subplots_adjust(left=0.14,bottom=0.15) #Adjust spacing to prevent clipping of x and y labels
#plt.savefig('Signal_Ratio.svg')
Explanation: <font size="4"><p>and plot the results:</p>
End of explanation
<END_TASK>