import openpyxl as px # Excel .xlsxファイルを扱うモジュールを読み込む
from math import exp
import numpy as np
from matplotlib import pyplot as plt
"""
Excel .xlsxファイルから太陽光スペクトルを読み込み、Planckの公式と比較する
"""
# 定数
pi = 3.14159265358979323846
h = 6.6260755e-34 # Js";
hbar = 1.05459e-34 # "Js";
c = 2.99792458e8 # m/s";
e = 1.60218e-19 # C";
kB = 1.380658e-23 # JK-1";
# 太陽光スペクトルの xlsxファイル
spectrum_file = 'solar.xlsx'
# 黒体の温度
T = 6000 # K
# 実行時引数を処理
argv = sys.argv
if len(argv) == 1:
print("Usage: python blackbody_radiation.py T")
print(" e.g.: pythnon blackbody_radiation.py 6000")
exit()
if len(argv) >= 2:
T = float(argv[1])
def main():
# Excelファイルを読み込む。
book = px.load_workbook(spectrum_file)
name = book.get_sheet_names()
print(name)
# アクティブシートの内容を now変数に読み込む
now = book.active
x = []
y = []
i = 0
print("lambda(um)\tE(W/m2/um)")
while 1:
# 読み込むcellの指定。B,D行はAM1.5 直達光、H,J行はAM1.5 全天光
# データは2行目から始まるので、初期値が i=0 の場合は i+2 で指定する
# xcell = "B" + str(i+2)
# ycell = "D" + str(i+2)
xcell = "H" + str(i+2)
ycell = "J" + str(i+2)
# cellの内容をfloat()で浮動小数点に変換。失敗したらデータ終了と判断
try:
x.append(float(now[xcell].value) * 1.0e3) # um => nm
y.append(float(now[ycell].value) * 1.0e-3) # W/m2 um => W/m2 nm
except:
break
print("%8.5g\t%8.5g" % (x[i], y[i]))
i += 1
# Planckの公式から計算
yp = []
for i in range(len(x)):
wl = x[i] * 1.0e-9
k = 8.0 * pi / pow(wl, 5) * h * c
expx = exp(h * c / kB / T / wl)
E = k / (expx - 1.0) * 1.0e-6 * 1.5 # W/(m2 m) => W/(m2 nm)
yp.append(E)
#=============================
# グラフの表示
#=============================
fig = plt.figure()
ax1 = fig.add_subplot(1, 1, 1)
ax1.plot(x, y, label = spectrum_file)
ax1.plot(x, yp, label = 'Planck eq', color = 'red')
ax1.set_xlabel("wavelength (um)")
ax1.set_ylabel("Energy density (W/m2nm)")
ax1.set_xlim([0, 3000.0])
ax1.set_ylim([0, 2.2])
ax1.legend()
plt.pause(0.1)
print("Press ENTER to exit>>", end = '')
input()
if __name__ == '__main__':
main()
