How to create an inverse pole figure color map?

Edit: I could simplify the code even more by accessing the 3D cooridnates directly, instead of using the inverse projection.

---

Disclaimer

As I finally understand what's going on, I thought I'd present my own solution and give some explanations. It's based on the dame idea as Rahul Narain's answer, only a bit simpler. The example image I provided in the question was a very poor choice, in that the correct result looks a little different (and doesn't need the parameters a and b, that Rahul introduced).

Anyway, I hope this will maybe help someone who is just as much at loss with (inverse) pole figures as I used to be a couple of weeks ago.

Basic idea

When dealing with cubic structures, every direction can be transformed to an equivalent one, whose stereographic projection lies within the stereographic standard triangle. As an example, the following picture shows the stereographic projection of some specific directions (integer points) onto the plane $z=-1$. The orange dots represent directions of the form $\langle 4,2,5 \rangle$, i.e. directions equivalent to $[4,2,5]$ by symmetry. 24 more such directions are projected to outside of the unit circle and thus are not visible in this plot. The green dot shows the equivalent direction whose projection lies inside the standard triangle.

The corners of the standard triangle are directions of the kind $\langle 1,0,0 \rangle$ (bottom left), $\langle 1,1,0 \rangle$ (bottom right) and $\langle 1,1,1 \rangle$ (top right), respectively. (In my case, to be more precise, $[0,0,1]$, $[0,1,1]$, and $[1,1,1]$.) Each direction $[x,y,z]$ can be expressed by a triple $(u,v,w)$ such that $$[x,y,z]= u \cdot[0,0,1]+v \cdot[0,1,1]+ w \cdot[1,1,1].$$ These are obviously: $$u=z - y,\quad v=y - x,\quad w=x .$$ This triple is used to assign each direction inside the standard triangle a unique RGBColor.

Step by step code

First, I define the stereographic projection:

project[v : {x_, y_, z_}] := {x, y}/(Norm@v + z)

The standard triangle can be easily parametrized by:

ParametricPlot[project@{s, t, 1}, {s, 0, 1}, {t, 0, s}]

The ColorFunction is similar to the one presented by Rahul and uses $u$, $v$, and $w$ from above:

makeColor[{x_, y_, z_}] := 
 RGBColor @@ ({z - y, y - x, x}/Max@{z - y, y - x, x})

This makeColor function has to be applied to the 3D coordinates, which can be accessed via the parameters of the ParametricPlot, i.e. #3 and #4:

ColorFunction -> (makeColor[{#4, #3, 1}] &)

Everything put together

Finally, everything is put together with some additional options for ParametricPlot:

project[v : {x_, y_, z_}] := {x, y}/(Norm@v + z)
makeColor[{x_, y_, z_}] := 
 RGBColor @@ ({z - y, y - x, x}/Max@{z - y, y - x, x})

ParametricPlot[project@{s, t, 1}, {s, 0, 1}, {t, 0, s}, 
 ColorFunction -> (makeColor[{#4, #3, 1}] &),
 Mesh -> None,
 BoundaryStyle -> Black,
 Frame -> False,
 Axes -> False]

I'm just going to leave this here. It started out as an attempt at a mathematically reasonable colouring of an orthographic projection of the fundamental domain, but that doesn't resemble the figure you want, so then I was forced to add a couple of artistic knobs.

z[x_, y_] := Sqrt[1 - x^2 - y^2]
normalize[RGBColor[r_, g_, b_]] := RGBColor @@ ({r, g, b}/Max[r, g, b])
a = 1.2;
b = 1.3;
ParametricPlot[{x, y}, {x, 0, 1/Sqrt@2}, {y, 0, 1/Sqrt@3}, 
 RegionFunction -> Function[{x, y}, y <= x && x <= z[x, y]], 
 Mesh -> None, Frame -> None, Axes -> None, 
 ColorFunction -> 
  Function[{x, y}, 
   normalize[
    RGBColor[1 - ArcTan[z[x, y], x]/(Pi/4), 
     a ArcTan[z[x, y], (x - y)/Sqrt@2]/ArcTan[1/Sqrt@2], 
     b ArcTan[z[x, y], y]/(Pi/4)]]], 
 ColorFunctionScaling -> False, PlotPoints -> 100]

You can tweak the values of a and b to adjust the position of the white point.

I do find the solution proposed by einbandi very elegant and very simple, but it places the white point at hkl coordinates <1,2,3> and this does not correspond to that used conventionally (see the image posted in the initial question).

Let us assume that the point of which we want the color is called O and has its stereographic projection as in the figure below.

Here is the rule used conventionally: The position of the point R, G, B is:

X(R)=0
Y(R)=0
X(G)=sqrt(2)-1
Y(G)=0
X(B)=sqrt(3)/2-1/2
Y(B)=X(B)

The intersection point E of a line AB with a line CD is given by:

 K1=(X(A)-X(B))*(Y(C)-Y(D))-(Y(A)-Y(B))*(X(C)-X(D))
 X(E)=((X(A)*Y(B)-Y(A)*X(B))*(X(C)-X(D))-(X(A)-X(B))*(X(C)*Y(D)-Y(C)*X(D)))/K1
 Y(E)=((X(A)*Y(B)-Y(A)*X(B))*(Y(C)-Y(D))-(Y(A)-Y(B))*(X(C)*Y(D)-Y(C)*X(D)))/K1

Then the position of G' is the intersection of GO and RB, and the position of B' is the intersection of BO and RG.

The position of R' cannot be calculated with the same procedure since BG is a portion of circle centered at [-1,0] and of radius sqrt(2). It is given by:

K1=(Y(O)/X(O))^2
X(R')=(sqrt(K1+2)-1)/(K1+1)
Y(R')=Y(O)/X(O)*X(R')

Finally, the level of each color is given by

• red = ||OR'||/||RR'||
• green = ||OG'||/||GG'||
• blue = ||OB'||/||BB'||

Then each color is normalized by max(red, blue, green)

PS: A came randomly on this page (Google) but I'm not a Mathematica user. So feel welcome to convert my method into Mathematica code.