# Need help simulating solar limb darkening

edit: At the Solar Dynamics Observatory (SDO) website, I just found the image sdo.gsfc.nasa.gov/assets/img/latest/latest_1024_HMIIC.jpg. The color gradient of the limb darkening seems very similar to the Wikimedia image below. I've discovered that it is called a "colorized intensitygram" and the color gradient is purely artificial - the data is single channel intensity. The limb darkening is certainly real (compare to the artificially "flattented" display!)

I would like to try to start to understand what the sun actually "looks like", with human vision, assuming the brightness has been reduced.

Like the gas giant planets, the sun doesn't have an abrupt surface, it just gets denser and hotter and denser and hotter as you go deeper.

One consequence of this is limb-darkening. As the material becomes denser the deeper you go, it becomes more opaque, so - roughly speaking - the light you see, including the color and brightness, is determined by the layers above that point.

If you look at center of the solar disk, you can see deeper and hotter parts. If you look near the edge of the sun, or solar limb, the incidence is oblique and you are not seeing as deep or as hot.

The solar physics and photon transport theory is complicated, but in the visible part of the spectrum it may be OK to think of blue light scattering much more strongly than red light. In the thick, dense solar atmosphere, it scatters enough to attenuate. So in blue light you are seeing even shallower, which is even colder, and therefore dimmer in blue.

Right now I'd just like some good analytical approximation of the wavelength dependent limb darkening of the sun, or some actual linear images (before web-processing) of the sun at various visible wavelengths so I can make my own approximation.

Here is an image from Wikimedia titled 2012_Transit_of_Venus_from_SF which I've separated into RGB channels and plotted brightness across horizontal (solid) and vertical (dashed) diameters, 20 pixels wide. You can see the dramatic difference in behavior of different wavelengths. Since the image is "from the internet" I've no information about the linearity, so it's an illustrative example, but not what you'd think of as data.

import numpy as np
import matplotlib.pyplot as plt

w, c, hw  = 600, 300, 10   # image happens to be 600x600 pixels square!

hor = img[    c-hw:c+hw].sum(axis=0) / (2*hw)
ver = img[ :, c-hw:c+hw].sum(axis=1) / (2*hw)

plt.figure()

plt.imshow(img)

r, g, b = c - 0.67*c*hor.T

plt.plot(r, '-r')
plt.plot(g, '-g')
plt.plot(b, '-b')

r, g, b = c - 0.67*c*ver.T

plt.plot(r, '--r')
plt.plot(g, '--g')
plt.plot(b, '--b')

plt.plot([0, w], [c, c], '-k')
plt.xlim(0, 599)
plt.ylim(599, 0)

plt.text(150, 137, 'R', fontsize=18)
plt.text(156, 192, 'G', fontsize=18)
plt.text(170, 272, 'B', fontsize=18)

plt.show()

• You are right, that image is bad. The true sun is white, so anything you see is off. Haven't found a good source for such an image, however... – PearsonArtPhoto Jun 8 '16 at 16:13
• @PearsonArtPhoto Over the years I've probably read Billmeyer and Saltzman's Principles of Color Technology cover-to-cover at least twice. While I understand that both "the sun is white" and "white is: the color of the sun" can be the right way to look at it, and brown really IS just dark yellow, I still think I'll have better luck figuring out manifolds and bifurcations before I'll ever be able to understand color spaces. – uhoh Jun 8 '16 at 18:38
• @PearsonArtPhoto I'm hoping there may be a satellite, or a telescopic study, with linear monochrome images in wavelength bands in the visible range. That, or an analytical approximation based on that kind of data in an old paper, written back when it would have been simple enough to read - without the magnetohydrodynamics etc.. – uhoh Jun 8 '16 at 18:42
• @PearsonArtPhoto but while you might call the color of sunlight white, you would then be forced to admit the middle of the sun was slightly blue to offset the reddening of the solar limb. There really is something like an effective color temperature shift from center to edge, – uhoh Jun 8 '16 at 18:45
• Possibly helpful: Wavelength dependency of the Solar limb darkening – 2012rcampion Jun 9 '16 at 15:50

I used the paper Wavelength dependency of the Solar limb darkening for solar limb darkening data.

It uses the following model for the normalized brightness distribution across the disk of the Sun:

$$I(\mu)=1-u(1-\mu^\alpha)$$

Here, $\mu$ is the normalized distance from the limb; when far from the Sun, it can be expressed in terms of the normalized distance from the center of the disk, $r$:

$$\mu = \sqrt{1-r^2}$$

$u$ and $\alpha$ are parameters. The paper is not very clear about the value of $u$, but it seems like they used $u=0.85$.

This formula tells us how the brightness varies from the center of the solar disk to the limb for a given wavelength. In order to properly represent the color of the Sun, we also need to know the relative brightness of different wavelength bands.

For this we need data on the solar spectrum, which I obtained from here (this data is for "air mass zero," meaning our result will be what the Sun would look like outside the Earth's atmosphere).

I assume that this spectrum is averaged over the whole disk, in which case we need to be able to scale the limb darkening model by average instead of peak (center) brightness. The formula given in the paper is:

$$I_\text{avg} \propto \int_0^1 I(\mu)\mu\ d\mu \propto \frac{2 + \alpha (1 - u)}{2 + \alpha}$$

The last piece of the puzzle is to convert the spectrum to a color. In this case I first compute the color in the CIE XYZ color space, then transform to linear RGB and then sRGB. (See this answer I posted on photography.SE for more details about computing the XYZ color of a spectrum, and the sRGB Wikipedia page for more detail about the conversion to sRGB.) I used the standard colorimetric observer tables from here.

Here is the Mathematica code that I used:

alphaRaw = ImportString[(* copy and paste from the PDF table *), "Table"];
(* join corresponding columns *)
Flatten[Partition[%, 2], {{2}, {1, 3}}];
(* resample to 1 nm *)
{alphaPS, alphaNL} = (Interpolation[#] /@ Range[360, 830]) & /@ %;


spectrumRaw = Import["http://rredc.nrel.gov/solar/spectra/am0/E490_00a_AM0.xls"];
(* resample to 1 nm (inputs in um!) *)
spectrum = Interpolation[spectrumRaw[[1, 2 ;;, ;; 2]]] /@ (Range[360, 830] / 1000);


cieRaw = Import["http://www.cis.rit.edu/research/mcsl2/online/CIE/StdObsFuncs.xls"];
(* pair each color function column with the wavelength column *)
Transpose[{cieRaw[[1, 6 ;;, 2 ;; 4]]], #}] & /@ cieRaw[[1, 6 ;;, 2 ;; 4]]];
(* resample to 1 nm *)
cieXYZ = (Interpolation[#] /@ Range[360, 830]) & /@ %;


With[{u = 0.85},
(* compute color at 101 points evenly spaced in mu *)
colorXYZ = Transpose @ Table[
Total[ (* integration -> summation *)
Transpose @ cieXYZ
* spectrum
* (1 - u (1 - mu^alphaPS))
/ ((2 + alphaPS (1 - u)) / (2 + alphaPS))
],
{mu, 0, 1, 0.01}
];
(* scale so maximum = 1 *)
colorXYZ /= Max[colorXYZ / {0.97, 1, 0.83}];
(* dim by 15% *)
colorXYZ *= 0.85;
]


(I'm not sure why Mathematica renders the XYZ colors with a blue tint.)

colorRGBlinear = {{3.2406, -1.5372, -0.4986}, {-0.9689, 1.8758, 0.0415}, {0.0557, -0.2040, 1.0570}}.colorXYZ;
colorRGBlinear /= Max[colorRGBlinear];
colorRGBlinear *= 0.85;


gamma = With[{a = 0.055}, If[# <= 0.0031308, 12.92 #, (1 + a) #^(1/2.4) - a] &];
colorSRGB = Map[gamma, colorRGBlinear, {2}];


Plot of the color difference (in CIE xy coordinates) from the limb to center (small black arrow). The ticks on the spectral locus are wavelength in nm, and the ticks on the Planckian locus are temperature in K.

• Nice!! You know those CIE color space diagrams with a Color Temperature line? Are you able to make one with dots at say 0, 0.5, 0.8, 0.9, and 0.95 $r_{disk}$ and add a "White" dot? Or some other way to show a slight color shift for people who may not have sufficient "hue acuity" to see it here? And/or make some kind of diameter R, G, and B plot? These images of "what the sun would look like" are beautiful for sure - but I need something quantitative. In the mean time I'll try to post a Python equivalent for those who don't have US $300 – uhoh Jun 10 '16 at 3:20 • @uhoh The colors stay very close to white, with a color temperature of around 5100 K at the limb and 6200 K at the center. The xyY plot isn't very interesting, it just looks like a tiny segment of the Planckian locus. As for the atmosphere, I don't think there will be any effect (the linear RGB picture looks pretty much identical to what I've seen in solar telescopes). – 2012rcampion Jun 10 '16 at 4:10 • OK could you find at least some way to show your center to edge shift that would still be accessible to the more than 4% of stackexchange users who have some form of color-blindness? I've suggested two methods, maybe there's another, but please, something unambiguous, in addition to the barely detectable hue gradient pattern. Thanks! note: I'm updating my plot in the question as well :) – uhoh Jun 10 '16 at 4:30 • @uhoh What is the end result you're trying to achieve? Would a table be more helpful than a plot? – 2012rcampion Jun 10 '16 at 4:39 • For me personally, I can follow your explanation, math, and Mathematica just fine - it's GREAT! (The paper also has some analytical expressions for the falloff of$\alpha\$ with wavelength as a short-cut). I think your answer is also wonderful for the general reader with "nominal" color perception. I'm just asking for one more thing: a non-hue-perception-based display of color shift vs radius that a stackexchange user with monochrome vision could still make use of. Like a plot or graph, "engineer style" :) – uhoh Jun 10 '16 at 4:55