プランクの法則を綺麗にプロットしたい

技術的に難しい部分とか特に無いのでこちらで。

プランクの法則

Wikipediaで確認しよう。ちなみに英語版のほうが詳しいぞ。

プランクの法則 - Wikipedia

あまりきちんと説明するとボロがでるのでいい加減に説明すると、暖かいものほど強い光を出して、暖かいものほど光の波長のピークが短波長側にずれるというやつだ。

よくあるプロット

こういうのよく見かけるだろう。

f:id:Rion778:20170502233943p:plain

この手のプロットだと英語版Wikipedia(Planck’s law - Wikipedia)にあるやつが割と好きだ。太陽の表面温度くらいの放射スペクトルがちょうど可視光になってるのがよく分かる。

ほしいプロット

温度も連続変数のはずだ。だとすれば、なんかこう、グラデーションっぽい感じでプロットしてほしい。

やっていく

パッケージ

今回はggplot2だけです。

library(ggplot2)

計算

k = 1.38e-23 # Boltzmann constant
h = 6.62e-34 # Planck constant
c = 3e8      # speed of light

p_rad <- function(lambda, t){
  2*h*c^2/(lambda^5) * 1/(exp((h*c)/(lambda*k*t)) - 1)
}

lim = 5000
t = lim - (1:(lim^(1/1.3)))^1.3
lambda = 10^-9 * 1:2000
param = expand.grid(t=t, lambda=lambda)
u = cbind(param, rad = p_rad(param$lambda, param$t))

温度と波長について少しずつ値を変えながら色々な組合せでの放射強度を計算していくわけだけれども、温度については先程のプロットからも分かるように、高温側では少し値が変わっただけで放射強度も大きく変わり、低温側ではそれほど変わらないという関係がある。

従って、高温側の隙間を小さく、低温側の隙間を大きくするような感じで調整してやると多分ポイントのムラが減る。

また、「色々な組み合わせ」はexpand.gridでやればパッと作れる。最初forとrbindで適当にやったら大変時間がかかった。

プロット

jet.colors <- colorRampPalette(c("#00007F", "blue", "#007FFF", "cyan", "#7FFF7F", "yellow", "#FF7F00", "red", "#7F0000"))

ggplot(u, aes(x = lambda, y = rad, color = t)) +
  geom_point(size = 0.5) + 
  scale_color_gradientn(name = "Temperature(K)", colors = jet.colors(7)) +
  theme_classic() +
  labs(x = "Wavelength (m)", y = "Intensity (W/sr/m^2)") +
  theme(legend.position = "right")