Я пытаюсь изобразить черное тело на разных реперах на одном графике?

Question

Я пытаюсь изобразить черное тело при разных значениях в одном и том же графике, температуры составляют 10, 10 ^ 3, 10 ^ 5, 10 ^ 7, 10 ^ 9 и 10 ^ 12 Кельвин. 10 ^ 12 выглядит хорошо, но остальные выглядят плоскими. Я уже пытался поместить его в масштаб журнала, но, похоже, это не помогло.

Вот мой код:

import numpy as np 
import matplotlib.pylab as plt

k = 1.38 * 10**-23  

h = 6.62607015e-34

c = 3e8

def pla(T, nu):
    x = (h*nu) / (k*T)

    c1 = 2 * h * nu**3

    planck = c1 / (np.exp(x) * c**2)

    return planck

nus = np.linspace(0, 4e23)

T10 = pla(10, nus)
Te3 = pla(1e3, nus)
Te5 = pla(1e5, nus)
Te7 = pla(1e7, nus)
Te9 = pla(1e9, nus)
Te12 = pla(1e12, nus)

plt.xlabel("frecuencia")
plt.ylabel("temperatura")

plt.plot(nus, T10, color='blue', label='T10')
plt.plot(nus, Te3, color='red', label='Te3')
plt.plot(nus, Te5, color='orange', label='Te5')
plt.plot(nus, Te7, color='purple', label='Te7')
plt.plot(nus, Te9, color='gray', label='Te9')
plt.plot(nus, Te12, color='green', label='Te12')

plt.legend(loc='upper right')

plt.show()

Результат:

Answer 1

Вот реализация, которую я нашел по адресу: { ссылка }

import matplotlib.pyplot as plt
import numpy as np

h = 6.626e-34
c = 3.0e+8
k = 1.38e-23

def planck(wav, T):
    a = 2.0*h*c**2
    b = h*c/(wav*k*T)
    intensity = a/ ( (wav**5) * (np.exp(b) - 1.0) )
    return intensity

# generate x-axis in increments from 1nm to 3 micrometer in 1 nm increments
# starting at 1 nm to avoid wav = 0, which would result in division by zero.
wavelengths = np.arange(1e-9, 3e-6, 1e-9) 

# intensity at 4000K, 5000K, 6000K, 7000K
intensity4000 = planck(wavelengths, 4000.)
intensity5000 = planck(wavelengths, 5000.)
intensity6000 = planck(wavelengths, 6000.)
intensity7000 = planck(wavelengths, 7000.)

plt.plot(wavelengths*1e9, intensity4000, 'r-') 
# plot intensity4000 versus wavelength in nm as a red line
plt.plot(wavelengths*1e9, intensity5000, 'g-') # 5000K green line
plt.plot(wavelengths*1e9, intensity6000, 'b-') # 6000K blue line
plt.plot(wavelengths*1e9, intensity7000, 'k-') # 7000K black line

# show the plot
plt.show()