1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
|
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Wed Nov 5 12:58:59 2025
@author: christian
"""
"""
Problem:
- Import xls data into Python
- Plot the Intensity [a.u.] vs pixels
- Interpolate and convert x-axis from pixels to nm (true wavelength) using Hg lamp data (using data in file: **Lampa_Calibrare_Mercur.xlsx**)
- Find response function of the spectrometer using the tungsten lamp data from file: "**Calibrare Intensitate Oxigen.xlsx**)": $R=\frac{I_{measured}}{I_{true}}$ (where True is computed by Planck's law of radiation (see notes in the pptx above)
- Convert y-axis from Intensity [a.u.] into Intensity in [W/(cm^2*sr*nm)] by dividing the measured Oxygen spectrum with the response function: $I_{oxygen, true}=\frac{I_{oxygen, measured}}{R}$
- Once the spectra is in real units: compute the density of one of the oxygen lines by integrating underneath one of the peaks (see equation from Slide 39 - bottom). We will give all of the constants that are in this equation (see the "Intensity_Calibration_Oxygen_Discharge_Solution.xlsx")
"""
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
Hg_data = 'Lampa_Calibrare_Mercur.xlsx'
#Ar_data = 'Spectru Descarcare Argon.xlsx'
Ox_data = 'Calibrare Intensitate Oxigen.xlsx'
df_Hg = pd.read_excel(Hg_data, header=11, engine='openpyxl')
df_Ox = pd.read_excel(Ox_data, header=3, engine='openpyxl')
#df_Ar = pd.read_excel(Ar_data, header=14, engine='openpyxl')
# Trims mercury data to the first 2000 pixels
df_Hg = df_Hg[df_Hg['Pixels']<2000]
# %% Plot Intensity vs Pixels for mercury
# Add Vertical lines to help match pixel count
wavelength_true = np.array([365.015, 404.656, 407.783, 435.833, 546.074, 576.960, 579.006])
bounds = np.array([635,1451])
m = (bounds[-1]-bounds[0])/(wavelength_true[-1]-wavelength_true[0])
c = bounds[0]-(m*wavelength_true[0])
pixel_count = m * wavelength_true + c
# Plot Hg intensity-pixel plot
plt.figure(figsize=(8,5))
plt.plot(df_Hg['Pixels'], df_Hg['Intensity'], linestyle='-')
plt.xlabel('Pixels')
plt.ylabel('Intensity [a.u.]')
plt.vlines(pixel_count,0,17000,colors='r')
plt.title('Pixel-Intensity (Hg)')
plt.grid(True)
plt.show()
# %% Calibrate length dimension using Mercury
# Arrays to be interpolated
wavelength_peak_nm = np.array([365.015, 404.656, 435.833, 546.074, 576.960, 579.006]) # From Literature
pixel_peak = np.array([635, 784, 902, 1324, 1443, 1451])
#Polynomial coefficients
C = np.polyfit(pixel_peak,wavelength_peak_nm,3)
# Calibration function
wavelength_Hg = lambda p : C[3] + C[2]*p+C[1]*p**2+C[0]*p**3
# or
#def wavelength_Hg(p):
# return C[3] + C[2]*p+C[1]*p**2+C[0]*p**3
wavelength_calibrated = wavelength_Hg(df_Hg['Pixels'])
# Plot Hg intensity-wavelength plot
plt.figure(figsize=(8,5))
plt.plot(wavelength_calibrated, df_Hg['Intensity'], linestyle='-')
plt.xlabel('Wavelength [nm]')
plt.ylabel('Intensity [a.u.]')
plt.title('Wavelength-Intensity (Hg)')
plt.grid(True)
plt.show()
# %% Calibrate intensity dimension using Tungsten
# Now that we have calibrated the wavelength dimension can can now use
# a tungsten source to calibrate intensity.
wavelength = df_Ox['Wavelength [nm]']
I_W_measured = df_Ox['I_Tungsten [a.u.]']-df_Ox['I_Background [a.u.]']
# Calculate theoretical Intensity of Tungsten using Planks eqn. (TO BE CHECKED)
epsilon = 0.4
h = 6.626e-34 # J/s
c = 2.99e8 # m/s
k = 1.381e-23 # J/K
T = 1800 # K
I_planck = lambda wavelength0 : epsilon*((2*h*c**2)/(wavelength0**5)*1/(np.exp(h*c/(k*T))-1))
I_W_true = I_planck(wavelength)
#Response Function
R=I_W_measured/I_W_true
#I_plasma_meas=R*I_plasma_true
# Plot Response function vs lambda
plt.figure(figsize=(8,5))
plt.plot(wavelength, R, linestyle='-')
plt.xlabel('Wavelength [nm]')
plt.ylabel('R')
plt.title('Response function')
plt.grid(True)
plt.show()
# Calculate Intensity of Oxygen
I_Ox_measured = df_Ox['I_Plasma Oxygen [a.u.]']-df_Ox['I_Background [a.u.]']
I_Ox_true = R * I_Ox_measured
# Plot Hg spectra [a.u.] plot
plt.figure(figsize=(8,5))
plt.plot(wavelength, I_Ox_measured, linestyle='-')
plt.xlabel('Wavelength [$nm$]')
plt.ylabel('Intensity [a.u,]')
plt.xlim(750,870)
plt.title('Wavelength-Intensity ($Ox_{measured}$)')
plt.grid(True)
plt.show()
# Plot Hg spectra real units plot
plt.figure(figsize=(8,5))
plt.plot(wavelength, I_Ox_true, linestyle='-')
plt.xlabel('Wavelength [$nm$]')
plt.ylabel('Intensity [$W/m^2/sr/nm$]')
plt.xlim(750,870)
plt.title('Wavelength-Intensity ($Ox_{true}$)')
plt.grid(True)
plt.show()
# %% Calculate Dencity of Ox at upper state
from scipy.integrate import simpson
I_Ox_true = I_Ox_true.to_numpy()
# Integrate numerator
delta = 1
peak = 846.5
bound_lower = peak - delta
bound_upper = peak + delta
mask = (wavelength >= bound_lower) & (wavelength <= bound_upper)
wavelength_trimmed = wavelength[mask]
I_Ox_trimmed = I_Ox_true[mask]
numerator = simpson(I_Ox_trimmed, wavelength_trimmed)
v = c/(peak*10**(-9))
A = 3.6e6 # 1/s
l = 0.015 # m
n = numerator / (1/(4*np.pi)*h*v*A*l)
|
