35
$\begingroup$

Since ColorData["VisibleSpectrum"] is wrong, I would like to have a more accurate function to use.
Can this information be extracted from Mathematica itself?

Improve this question
$\endgroup$

3 Answers 3

Reset to default
36
$\begingroup$

Notice:

Simon Woods did just this months ago for an answer I missed:

It seems that it can. By spelunking ChromaticityPlot I found:

Image`ColorOperationsDump`$wavelengths Image`ColorOperationsDump`tris

These are a list of wavelengths and their corresponding XYZ color values used by this plot command:

ChromaticityPlot["sRGB", Appearance -> {"VisibleSpectrum", "Wavelengths" -> True}]

enter image description here

We can therefore use them to generate a new color function:

ChromaticityPlot; (* pre-load internals *) newVisibleSpectrum = With[ {colors = {Image`ColorOperationsDump`$wavelengths, XYZColor @@@ Image`ColorOperationsDump`tris}\[Transpose]}, Blend[colors, #] & ];

A comparison with the old function:

ArrayPlot[ {Range[385, 745]}, ImageSize -> 400, AspectRatio -> 0.1, ColorFunctionScaling -> False, ColorFunction -> # ] & /@ {"VisibleSpectrum", newVisibleSpectrum} // Column

enter image description here

589nm is now the bright sodium yellow that it should be:

Graphics[{newVisibleSpectrum @ 589, Disk[]}]

enter image description here

If you wish to integrate this into ColorData see:

As requested by J.M. red-green-blue plots for each function:

old = List @@@ Array[ColorData["VisibleSpectrum"], 361, 385]; new = List @@@ ColorConvert[Array[newVisibleSpectrum, 361, 385], "RGB"]; ListLinePlot[Transpose @ #, PlotStyle -> {Red, Green, Blue}, DataRange -> {385, 745} ] & /@ {old, new}

enter image description here

Clipping occurs during conversion to screen RGB; the newVisibleSpectrum function actually produces unclipped XYZColor data. For example projected over gray:

newVSgray = With[{colors = Array[{#, Blend[{newVisibleSpectrum@#, ColorConvert[GrayLevel[.66], "XYZ"]}, 0.715]} &, 361, 385]}, Blend[colors, #] &]; ListLinePlot[ List @@@ ColorConvert[Array[newVSgray, 361, 385], "RGB"] // Transpose, PlotStyle -> {Red, Green, Blue}, DataRange -> {385, 745}, ImageSize -> 400]

enter image description here

Which corresponds to the plot:

ArrayPlot[{Range[385, 745]}, ImageSize -> 400, AspectRatio -> 0.1, ColorFunctionScaling -> False, ColorFunction -> newVSgray, Background -> GrayLevel[0.567]]

enter image description here

cf. "VisibleSpectrum" similarly over gray blended in XYZColor and RGBColor respectively:

enter image description here
enter image description here

Note that neither rendering of this spectrum has the vibrancy of newVisibleSpectrum.

Improve this answer
$\endgroup$
17
  • 2
    $\begingroup$ It seems the internal functions are only available after you call ChromaticityPlot once. Is that normal? $\endgroup$ Commented Feb 3, 2015 at 13:40
  • $\begingroup$ @Rahul Yes, it's common with many internal functions and I should have thought to mention it. $\endgroup$ Commented Feb 3, 2015 at 14:12
  • 1
    $\begingroup$ +1 Nice. Though spelunking ChromaticityPlot is so last year :-) $\endgroup$ Commented Feb 3, 2015 at 14:27
  • $\begingroup$ @Simon Gah! How'd I miss that? $\endgroup$ Commented Feb 3, 2015 at 14:37
  • 1
    $\begingroup$ That gradient is very bad, and that makes me dubious about the black body one as well. I'll try to see where/how they are implemented and fix them. I should also read here more often... (btw if the color matching functions are such coveted material I may put them in a decent context or push for them to be in System altogether) $\endgroup$ Commented Aug 2, 2015 at 9:19
9
$\begingroup$

(with many thanks to halirutan and kirma for their kind assistance)

Here's a different take. In this article, piecewise Gaussian functions that approximate the CIE color matching functions are presented. For this answer, instead of just taking the coefficients from the paper directly, I used their proposed model in FindFit[], using a 1 nm tabulation of the 1931 CMFs from here as the data the proposed model is being fitted to. Here are the resulting functions:

SetAttributes[{cieX, cieY, cieZ}, Listable]; cieX[λ_] := {0.36263412, 1.05593554, -0.2122116}.Exp[ -MapThread[(λ - #1) Piecewise[{{#2, λ < #1}, {#3, True}}] &, {{440.76797, 599.44051, 493.87482}, {0.066588512, 0.027633445, 0.048888466}, {0.020005468, 0.032068622, 0.039095442}}]^2/2] cieY[λ_] := {0.82087906, 0.28579841}.Exp[ -MapThread[(λ - #1) Piecewise[{{#2, λ < #1}, {#3, True}}] &, {{568.78966, 530.87379}, {0.021328878, 0.061294726}, {0.024693511, 0.032178136}}]^2/2] cieZ[λ_] := {1.21649968, 0.68126275}.Exp[ -MapThread[(λ - #1) Piecewise[{{#2, λ < #1}, {#3, True}}] &, {{436.96247, 459.03433}, {0.084453906, 0.038519512}, {0.027788061, 0.072502761}}]^2/2]

In version 10, one can of course directly do Through[XYZColor[cieX, cieY, cieZ][λ]]; for older versions, you will have to do the conversion to sRGB yourself, like so:

(* gamma correction *) sRGBGamma = Function[x, With[{z = Abs[x]}, Sign[x] Piecewise[{{12.92 z, z <= 0.0031308}}, 1.055 z^(1/2.4) - 0.055]], Listable]; myVisibleSpectrum[λ_] := RGBColor @@ Clip[sRGBGamma[{{3.2404542, -1.5371385, -0.49853141}, {-0.96926603, 1.8760108, 0.041556017}, {0.055643431, -0.20402591, 1.0572252}} . Through[{cieX, cieY, cieZ}[λ]]], {0, 1}]

Here are the corresponding gradient and component plots for myVisibleSpectrum[]: myVisibleSpectrum demonstration

Here's the result of coloring a Disk[] with myVisibleSpectrum[589]: sodium yellow?

The approximation looks good, at least with my eyes.

This replicates the internal data (as used in newVisibleSpectrum) quite well:

cieXYZ = Through[XYZColor[cieX, cieY, cieZ][#]] &; {d1, d2} = Transpose /@ Table[List @@ fn[λ], {fn, {cieXYZ, newVisibleSpectrum}}, {λ, 385, 745}]; Show[ListLinePlot[d1, PlotStyle -> Directive[Thickness[0.01], Hue[0.55, 0.5, 1]], PlotTheme -> "Pastel", InterpolationOrder -> 3], ListLinePlot[d2, PlotStyle -> Directive[Black, Dashed]]]

comparison of tabulated CMF with analytic fit

Alternatives

I had already discussed Bruton's algorithm in my other answer. I'm not sure how to convert wavelengths to the colors used in Sam's answer, but here is how to see them in Mathematica:

DensityPlot[x, {x, 0, 1}, {y, 0, 1/8}, AspectRatio -> Automatic, ColorFunction -> Function[t, RGBColor @@ Haversine[2 π t - {π, 5 π/3, π/3}]]]

Sam's color gradient

Improve this answer
$\endgroup$
2
  • $\begingroup$ Thanks for answering. It seems you've replicated the internal data quite well. I'll assume your gedanken mma isn't up to direct .PNG output yet so I included a plot to show this. $\endgroup$ Commented Jun 27, 2015 at 2:56
  • $\begingroup$ No, I'm using a borrowed computer for now. :) (though it only has version 8). Thanks for putting in the plot! $\endgroup$ Commented Jun 27, 2015 at 3:14
5
$\begingroup$

Here's my mathematical rainbow spectrum plot.

I'm not sure how accurate it is, and I did not attempt to calibrate it with measured data; but it's the smoothest, best looking spectrum I've seen, so I thought I would post it here.

rainbow spectrum

Large images: http://sam.nipl.net/pix/spectrum.png http://sam.nipl.net/pix/spectrum2.png

Reasoning:

  • All colours in the rainbow spectrum should be shown with equal intensity.
  • In RGB, yellow(1,1,0) has greater intensity than red(1,0,0) or green(0,1,0).
  • So, we need a less intense yellow for our rainbow. I used cosines.
  • The rainbow is an octave of light, from 400THz-800Thz, or 750nm-375nm.
  • We can see 12 rainbow colours, like the 12 musical notes.
  • The refractive index varies with the wavelength of light.
  • The rainbow is a uniform plot over wavelength (not freq. or "note").
  • Several "rainbow echoes" are seen, due to resonance with higher octaves.

Here is a picture of rainbow echoes, these are due I suppose to harmonic resonance:

rainbow with echoes

Here is the code, it's not in Mathematica (sorry!):

#!/usr/local/bin/cz -- use b Main() space(1900,100) origin(w_2, h_2) rainbow_init() for(x, 0.0, w) rainbow_colour(x/w) line(x, 0, x, h) paint() if args && cstr_eq(arg[0], "-s"): plot_12_tones() num rb_red_angle, rb_green_angle, rb_blue_angle num rb_red_power, rb_green_power, rb_blue_power rainbow_init() rb_red_angle = deg2rad(0) rb_green_angle = deg2rad(120) rb_blue_angle = deg2rad(-120) rb_red_power = 1 rb_green_power = 1 rb_blue_power = 1 colour rainbow_colour(num h) num y = pow(2, h)-1 num a = 2*pi*y num r = rb_red_power * (cos((a-rb_red_angle))+1)/2 num g = rb_green_power * (cos((a-rb_green_angle))+1)/2 num b = rb_blue_power * (cos((a-rb_blue_angle))+1)/2 return rgb(r, g, b) plot_12_tones() black() for(i, 0, 12+1): num freq = pow(2, i/12.0) num wl = 2 - 1/freq * 2 num x = wl*w x = iclamp(x, 0, w-1) line(x, 0, x, 2) line(x, h, x, h-3)
Improve this answer
$\endgroup$
5
  • 2
    $\begingroup$ Sam, welcome to Mathematica Stack Exchange. I appreciate you posting this despite what I am about to say: if you do not have Mathematica an/or are unfamiliar with its language it will be difficult for you to meaningfully participate here. If you choose to participate you should focus on clear English descriptions of algorithms or approaches rather than code in another language. $\endgroup$ Commented Feb 5, 2015 at 17:29
  • 2
    $\begingroup$ Subjectively the yellow in your photo looks more intense than the one in your smooth spectrum; I take it your intent is not to replicate the photo but instead to produce a pleasingly smooth gradient, correct? $\endgroup$ Commented Feb 5, 2015 at 17:30
  • 1
    $\begingroup$ Mr.Wizard, yes that's correct. Perhaps I'll have a go at translating it into Mathematica... $\endgroup$ Commented Feb 5, 2015 at 17:43
  • 3
    $\begingroup$ The "echoes" are called supernumerary bows. en.wikipedia.org/wiki/Rainbow#Supernumerary_rainbow $\endgroup$ Commented Jun 17, 2015 at 13:56
  • $\begingroup$ These supernumerary bows have nothing whatsoever to do with "harmonic resonance," at least as such words mean. $\endgroup$ Commented Apr 22, 2020 at 1:44

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.