Area-preserving projections from hexagonal and triangular domains to the sphere and applications to electron back-scatter diffraction pattern simulations

In the materials community, uniform sampling on a sphere (either $\mathbb S^2$ or $\mathbb S^3$) is used in a variety of applications. Polycrystalline materials are commonly employed in engineering applications; they consist of large numbers of grains, separated from each other by grain boundaries. An individual grain is characterized by the orientation of its crystal (Bravais) lattice with respect to an external reference frame, and this orientation is usually expressed as an Euler angle triplet (φ₁, Φ, φ₂) [1], as a unit quaternion, q [2], or by Rodrigues parametrization [3]. In numerical simulations, one often desires to generate random sets of grain orientations by means of a uniform sampling of Euler space, or of the rotation space SO(3) (which can be mapped onto $\mathbb S^3)$. While several methods of uniform sampling of SO(3) have been described in the literature (e.g. [4, 5]), there is still a need for a fast and efficient algorithm that will generate a uniform grid on $\mathbb S^3$, suitable for interpolations and for representation of functions expanded in generalized spherical harmonics.

Crystal orientations are usually extracted from electron back-scatter diffraction (EBSD) patterns, which can be obtained in a scanning electron microscope equipped with a special EBSD detector. A sub-micron diameter high energy electron beam is scanned across the sample, and dwells at individual locations for a preset amount of time; the sample is typically oriented at a 70° angle with respect to the incident beam, to maximize the number of back-scattered electrons that are intercepted by the area detector. The distribution of these electrons is not uniform, and is spatially modulated because of the interactions of the outgoing electrons with the crystal lattice; the resulting diffraction pattern on the detector shows sets of parallel lines (Kikuchi bands) of varying widths and orientations [6]. Analysis of the locations and orientations of these bands (usually by means of a Hough transform) then allows for the determination of the orientation of the crystal lattice at the illuminated point. When the electron beam is scanned across a sample region, repeated observation and analysis of the EBSD patterns then results in an orientation map, i.e., a map that shows the variations in crystal orientation as a function of position on the sample surface.

For the numerical simulation of EBSD patterns, one can pre-compute a 'master EBSD pattern' by considering a uniform grid of outgoing beam directions represented by direction cosines on $\mathbb S^2$, and computing the back-scatter yield for each orientation; this latter computation involves integrating the Schrödinger equation for an electron traveling through a periodic lattice potential [7]. An individual EBSD pattern can then be computed by interpolating from the grid on $\mathbb S^2$, or, more conveniently, by interpolating on a square 2-D grid obtained from the $\mathbb S^2$ grid by an area-preserving projection [8]. Amongst the seven traditional crystal systems (cubic, tetragonal, orthorhombic, hexagonal, trigonal, monoclinic and anorthic), all but the hexagonal and trigonal symmetries can be represented properly on a square grid. For the hexagonal and trigonal (or rhombohedral) crystal systems, EBSD pattern simulation would benefit from the availability of an area-preserving projection between the sphere $\mathbb S^2$ and a hexagonal or triangular grid instead of a square grid, since such non-square grids are fully compatible with the rotational and mirror symmetry elements of these two crystal systems. The present contribution derives such area-preserving projections and illustrates how they can be used in the context of dynamical EBSD pattern simulations.

2.1. Derivation of the area-preserving projections

Consider a regular hexagon H_a, and an equilateral triangle, T_b, with edge lengths a and b, respectively, and centered at the origin O, as illustrated in figure 1. For each shape, let $\mathcal C_R$ be a circle with radius R $(R_{H_a}$ and $R_{T_b})$ with the same area as H_a and T_b. Since the areas of the hexagon and the triangle are given by

$$\begin{equation} \mathcal A (H_a) = \case{3}{2}\sqrt{3}\,a^2\tqs {\rm and} \tqs \mathcal A (T_b) = \case{\sqrt{3}}{4}\,b^2, \end{equation} \tag{ 1 }$$

we have

$$\begin{equation} \fl R_H=\sqrt{\frac{3\sqrt{3}}{2\pi}}a \approx 0.909\,312\, a\tqs {\rm and}\tqs R_T=\sqrt{\frac{\sqrt{3}} {4\pi}}b\approx 0.371\,258\, b. \end{equation} \tag{ 2 }$$

Each shape is subdivided into triangular sectors. For the hexagon we consider the sextants I_k, k = 0, ..., 5, which are defined as follows (see figure 1):

$$\begin{eqnarray} I_0 &= \left\{(x,y)\in \mathbb{R}^2, 0\le x,\ -x/\sqrt{3}\le y\le x/\sqrt{3}\right\}, \end{eqnarray} \tag{ 3a }$$

$$\begin{eqnarray} I_1 &= \left\{(x,y)\in \mathbb{R}^2, 0\le x,\ x/\sqrt{3}\le y \right\}, \end{eqnarray} \tag{ 3b }$$

$$\begin{eqnarray} I_2 &= \left\{(x,y)\in \mathbb{R}^2, x\le 0,-x/\sqrt{3}\le y \right\}, \end{eqnarray} \tag{ 3c }$$

$$\begin{eqnarray} I_3 &= \left\{(x,y)\in \mathbb{R}^2, x\le 0,\ x/\sqrt{3}\le y\le -x\sqrt{3}\right\}, \end{eqnarray} \tag{ 3d }$$

$$\begin{eqnarray} I_4 &= \left\{(x,y)\in \mathbb{R}^2, x\le 0,\ y\le x/\sqrt{3} \right\}, \end{eqnarray} \tag{ 3e }$$

$$\begin{eqnarray} I_5 &= \left\{(x,y)\in \mathbb{R}^2, 0\le x,\ y\le -x/\sqrt{3} \right\}. \end{eqnarray} \tag{ 3f }$$

For (x, y) ∈ H_a ∩ I_k one has x ⩽ α if k = 0, 1, 5 and x ⩾ −α if k = 2, 3, 4, where $\alpha=a\sqrt 3/2$. For the equilateral triangle we consider the triangular sectors J_ℓ, ℓ = 0, 1, 2, defined by

$$\begin{eqnarray} J_0 &= \left\{(x,y)\in \mathbb{R}^2,\ -x\sqrt{3}\le y\le x\sqrt{3}\right\}, \end{eqnarray} \tag{ 4a }$$

$$\begin{eqnarray} J_1 &= \left\{(x,y)\in \mathbb{R}^2, \ 0\le y,\ x\sqrt{3}\le y\right\}, \end{eqnarray} \tag{ 4b }$$

$$\begin{eqnarray} J_2 &= \left\{(x,y)\in \mathbb{R}^2,\ y\le 0, \ y\le -x\sqrt{3} \right\}. \end{eqnarray} \tag{ 4c }$$

For (x, y) ∈ T_b ∩ J_ℓ one has x ⩽ β if ℓ = 0 and −2β ⩽ x ⩽ β if ℓ = 1, 2, where $\beta =b/(2\sqrt 3)$.

Zoom In Zoom Out Reset image size

In the following, we will define two functions $\mathcal T:\mathbb R^2 \to \mathbb R^2$, which map regular hexagons $(\mathcal T_H)$ and equilateral triangles $(\mathcal T_{T})$ onto circles and preserve areas, i.e.

$$\begin{equation} \mathcal A(D)=\mathcal A(\mathcal T(D)), \end{equation} \tag{ 5 }$$

for every $D\subseteq \mathbb R^2$.

In order to deduce the formulas for $\mathcal T$ we restrict the discussion for the moment to the dodecant I in figure 1(a), bounded by the lines of equations y = 0 and $y=x/\sqrt{3}$, and the upper-half J of the sector J₀ in figure 1(b), bounded by the lines of equations y = 0 and $y=x\sqrt{3}$.

Every half-line in between these bounding lines has the equation y = mx, with $0\le m\le 1/\sqrt{3}$ for the hexagon and $0\le m\le \sqrt{3}$ for the equilateral triangle. The mappings $\mathcal T_H$ and $\mathcal T_{T}$ will be defined in such a way that the half-line y = mx is mapped onto the half-line of equation

$$\begin{equation} y = \varphi(m)x. \end{equation} \tag{ 6 }$$

For the regular hexagon, φ satisfies the conditions

$$\begin{equation} \fl \varphi(0) =0,\quad \varphi\left(\frac{1}{\sqrt{3}}\right)=\frac{1}{\sqrt{3}},\quad {\rm and}\quad 0\le\varphi(m)\le\frac{1}{\sqrt{3}}\quad {\rm for}\ 0\le m\le \frac{1}{\sqrt{3}}, \end{equation} \tag{ 7a }$$

while for the equilateral triangle φ satisfies

$$\begin{equation} \fl \varphi(0) =0,\quad \varphi(\sqrt{3})=\sqrt{3},\quad{\rm and}\quad 0\le\varphi(m)\le\sqrt{3}\quad \mbox{for}\ 0\le m\le \sqrt{3}. \end{equation} \tag{ 7b }$$

These definitions ensure that the bounding half-lines are invariant under the mapping $\mathcal T$, and that all half-lines y = φ(m)x lie between the bounding half-lines.

Next, we want to ensure that the half-edges of the hexagon H_a and the triangle T_b are each mapped onto the circle arcs inside the respective half-sectors I and J. We take the point M with coordinates (x_M, y_M) = (x_M, m x_M) ∈ I. The line OM then intersects the hexagon's half edge at the point Q with coordinates (α, mα), where m = y_M/x_M (see figure 1(a)); for the equilateral triangle, the point Q has the coordinates (β, mβ) (figure 1(b)).

The mapping $\mathcal T$ (either $\mathcal T_H$ or $\mathcal T_{T})$ will map the point Q onto a point P with coordinates (x_P, φ(m)x_P), situated on the circle $\mathcal C_R$ $(\mathcal C_{R_H}$ or $\mathcal C_{R_T})$, so we must have

$$\begin{equation} x_P^2+\varphi^2(m) x_P^2 = R^2. \end{equation} \tag{ 8 }$$

Hence,

$$\begin{equation} (x_P,y_P) = \left(\frac{R}{\sqrt{1+\varphi^2(m)}} ,\frac{R\,\varphi(m)}{\sqrt{1+\varphi^2(m)}}\right). \end{equation} \tag{ 9 }$$

Next we denote $N=\mathcal T(M)$, and we choose N such that

$$\begin{equation} \frac{ON}{OP} = \frac{OM}{OQ}. \end{equation} \tag{ 10 }$$

A simple calculation shows that

$$\begin{eqnarray} ON &= x_N\sqrt{1+\varphi^2(m)}, \end{eqnarray} \tag{ 11a }$$

$$\begin{eqnarray} OM &= x_M\sqrt{1+m^2}, \end{eqnarray} \tag{ 11b }$$

$$\begin{eqnarray} OQ &= \gamma\sqrt{1+m^2}, \end{eqnarray} \tag{ 11c }$$

$$\begin{eqnarray} OP &= R, \end{eqnarray} \tag{ 11d }$$

with γ = α for the hexagon, and γ = β for the equilateral triangle. From condition (10) one can easily deduce that

$$\begin{equation} (x_N,y_N) = \left(\frac{x_M R}{\gamma\sqrt{1+\varphi^2(m)}},\frac{x_M R\,\varphi(m)}{\gamma\sqrt{1+\varphi^2(m)}}\right). \end{equation} \tag{ 12 }$$

In conclusion, the transformation $\mathcal T$ will map the point (x, y) ∈ I into (X, Y) ∈ I, given by

$$\begin{equation} (X,Y) = \left(\frac{x R}{\gamma\sqrt{1+\varphi^2(\frac{y}{x})}},\frac{x R\,\varphi(\frac{y}{x})}{\gamma\sqrt{1+\varphi^2(\frac{y}{x})}}\right). \end{equation} \tag{ 13 }$$

The function φ will be determined in such a way that $|{\bit J}(\mathcal T)|=1$, where the Jacobian ${\bit J}(\mathcal T)$ is defined as

$$\begin{equation} {\bit J}(\mathcal T):= \left | \begin{array}{@{}cc@{}} \frac{\partial X}{\partial x} & \frac{\partial X}{\partial y} \\ \frac{\partial Y}{\partial x} & \frac{\partial Y}{\partial y} \\ \end{array} \right |. \end{equation} \tag{ 14 }$$

A simple calculation shows that the condition |J(T)| = 1, which is necessary for an equal-area mapping, reduces to

$$\begin{equation} \frac{R^2}{\gamma^2}\frac{\varphi'(m)}{1+\varphi^2(m)}=1. \end{equation} \tag{ 15 }$$

The integration, together with the condition φ(0) = 0, gives

$$\begin{equation} \varphi(m) = \tan\frac{m\gamma^2}{R^2}. \end{equation} \tag{ 16 }$$

For the hexagon we have γ = α and $R=a(3\sqrt{3}/2\pi)^{1/2}$, so that

$$\begin{equation} \varphi_H(m) = \tan\frac{m\pi}{2\sqrt{3}}, \end{equation} \tag{ 17 }$$

and we can easily verify that

$$\begin{equation} \varphi_H\left(\frac{1}{\sqrt{3}}\right) = \tan\frac{\pi}{6}=\frac{1}{\sqrt{3}}, \end{equation} \tag{ 18 }$$

as expected. For the equilateral triangle, we have γ = β and $R=b(\sqrt{3}/4\pi)^{1/2}$, so that

$$\begin{equation} \varphi_{T}(m) = \tan\frac{m\pi}{3\sqrt{3}} \end{equation} \tag{ 19 }$$

and

$$\begin{equation} \varphi_{T}\Big(\sqrt{3}\Big) = \tan\frac{\pi}{3}=\sqrt{3}, \end{equation} \tag{ 20 }$$

as required.

Substituting φ(m) into (13), we finally obtain that the point (x, y) ∈ I (and respectively J) will be mapped into

$$\begin{equation} (X,Y)_H = 3^{\frac{1}{4}}\sqrt{\frac{2}{\pi}}x\left( \cos\frac{y\pi}{2\sqrt{3}x},\sin\frac{y\pi}{2\sqrt{3}x}\right) \end{equation} \tag{ 21a }$$

for the hexagon, and

$$\begin{equation} (X,Y)_{T} = \frac{3^{\frac{3}{4}}}{\sqrt{\pi}}x\left( \cos\frac{y\pi}{3\sqrt{3}x},\sin\frac{y\pi}{3\sqrt{3}x}\right) \end{equation} \tag{ 21b }$$

for the equilateral triangle. It is easy to prove that these relations are valid for the whole sectors I₀ and J₀, respectively. Next, we wish to extend the mappings $\mathcal T_H$ and $\mathcal T_{T}$ to the whole plane.

2.2. The hexagonal mapping $\mathcal{\rm T}_{\rm H}$

In this section, we will denote the restriction of $\mathcal T_H$ to the sextant I₀ by $\mathcal T_0$. Then we consider the rotation matrices

$$\begin{equation} \fl \mathcal R:= \left(\begin{array}{@{}cr@{}} \cos\frac{\pi}{3} & -\sin\frac{\pi}{3}\\ \sin\frac{\pi}{3} & \cos\frac{\pi}{3}\end{array}\right) \tqs {\rm and}\tqs \mathcal R_k:=\mathcal R^k=\left(\begin{array}{@{}cr@{}} \cos\frac{k\pi}{3} & -\sin\frac{k\pi}{3}\\ \sin\frac{k\pi}{3} & \cos\frac{k\pi}{3}\end{array}\right). \end{equation} \tag{ 22 }$$

One can see that, if M ∈ I_k, then $\mathcal R_k^{-1}(M)=\mathcal R_k^T (M) \in I_0$, and therefore the transformation $\mathcal T_H$ at the point M ∈ I_k can be defined as

$$\begin{equation} \mathcal T_H(M):= (\mathcal R_k \circ \mathcal T_0\circ \mathcal R_k^T)(M). \end{equation} \tag{ 23 }$$

Further, for (x, y) ∈ I₁ we have

$$\begin{eqnarray*}&&\mathcal R^T \cdot \left(\begin{array}{@{}c@{}}x\\ y\end{array}\right) &=&\frac 12 \left(\begin{array}{@{}c@{}} x+y\sqrt 3 \\ y-\sqrt 3x\end{array}\right)\in I_0,\\ \fl \mathcal T_0\Big( \frac{x+y\sqrt 3}2,\frac {y-\sqrt 3x}2\Big)&=&\frac{3^{\frac 14}}{\sqrt{2\pi}} (x+\sqrt 3y)\left( \cos \frac{\pi(y-\sqrt 3 x)}{2\sqrt 3(x+\sqrt 3 y)},\sin \frac{\pi(y-\sqrt 3 x)}{2\sqrt 3(x+\sqrt 3 y)} \right)\\ &=& \frac{3^{\frac 14}}{\sqrt{2\pi}} (x+\sqrt 3y) (\cos \omega,\sin \omega), \end{eqnarray*}$$

$$\begin{eqnarray*}\fl (\mathcal R \circ \mathcal T_0\circ \mathcal R^T)(x,y)&= &\frac{3^{\frac 14}}{\sqrt{2\pi}} (x+\sqrt 3y) \left(\cos\Big(\omega+\frac \pi 3\Big),\ \sin \Big(\omega+\frac \pi 3\Big)\right)\\ \fl &=& \frac{3^{\frac 14}}{\sqrt{2\pi}} (x+\sqrt 3y) \left( \sin \Big( \frac \pi 6-\omega \Big),\ \cos \Big( \frac \pi 6-\omega \Big) \right)\\ \fl &=& \frac{3^{\frac 14}}{\sqrt{2\pi}} (x+\sqrt 3y) \left( \sin\frac{2\pi x}{3(x+\sqrt{3}y)},\cos \frac{2\pi x}{3(x+\sqrt{3}y)}\right). \end{eqnarray*}$$

In the same manner as before, for (x, y) ∈ I₂ we get

$$\begin{equation} \fl (\mathcal R_2 \circ \mathcal T_0 \circ \mathcal R_2^T)(x,y)=\frac{3^{\frac 14}}{\sqrt{2\pi}} (x-\sqrt 3y) \left( \sin\frac{2\pi x}{3(x-\sqrt{3}y)},\ -\cos\frac{2\pi x}{3(x-\sqrt{3}y)}\right). \end{equation} \tag{ 24 }$$

Further, for the other sextants, we use the equalities

$$\begin{equation*} \mathcal R_3=-\mathcal I_2,\tqs \mathcal R_4=-\mathcal R_,\tqs \mathcal R_5=-\mathcal R_2=\mathcal R^T, \end{equation*}$$

so finally we obtain that $\mathcal T_H$ maps the point (x, y) ≠ (0, 0) into (X, Y)_H given by

• For (x, y) ∈ I₀ ∪ I₃,
  $$\begin{equation} \fl (X,Y)_H=3^{\frac 14}{\sqrt{\frac 2\pi}}\, x \left( \cos\frac{y\pi}{2\sqrt{3}x},\sin\frac{y\pi}{2\sqrt{3}x}\right); \end{equation} \tag{ 25a }$$
• For (x, y) ∈ I₁ ∪ I₄,
  $$\begin{equation} \fl (X,Y)_H = \frac{3^{\frac 14}}{\sqrt{2\pi}} (x+\sqrt 3y) \left( \sin\frac{2\pi x}{3(x+\sqrt{3}y)},\ \cos\frac{2\pi x}{3(x+\sqrt{3}y)}\right); \end{equation} \tag{ 25b }$$
• For (x, y) ∈ I₂ ∪ I₅,
  $$\begin{equation} \fl (X,Y)_H = \frac{3^{\frac 14}}{\sqrt{2\pi}} (x-\sqrt 3y) \left( \sin\frac{2\pi x}{3(x-\sqrt{3}y)},-\cos\frac{2\pi x}{3(x-\sqrt{3}y)}\right). \end{equation} \tag{ 25c }$$

For the origin we define $\mathcal T_H(0,0)=(0,0)$. Figures 2(a) and (b) show how a uniform hexagonal sampling grid inside a regular hexagon is converted to a uniform grid inside the unit circle. Note that grid points that lie on the edge of a hexagon inside the external hexagon (outlined in a thicker line) in (a) are mapped onto a single circle in (b). Figure 3 (left) shows a grid on a disc, obtained by dividing each triangle of the hexagon H_a into 36 small triangles of equal area.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

2.3. The equilateral triangle mapping $\mathcal{\rm T}_{{\rm T}}$

2.4. The inverse hexagonal map

To start, we restrict the discussion again to the first sextant I₀ and we determine the inverse of $\mathcal T_0$ given in (21a). In order to simplify the formulas, we denote $\theta={y\pi}/{2\sqrt 3 x}$ and $\eta=3^{\frac 14}{\sqrt{2/\pi}}$. With these notations, we have X² + Y² = η²x² and (Y/X) = tan θ, whence

$$\begin{equation} x=\frac{\sqrt{X^2+Y^2}}\eta \tqs {\rm and}\tqs y=\frac{2\sqrt 3}{\pi}\,x \arctan \frac YX. \end{equation} \tag{ 27 }$$

Therefore, $\mathcal T_0^{-1}$ maps a point (X, Y) ∈ I₀ into the point (x, y) ∈ I₀ given by

$$\begin{eqnarray} x&=& 3^{-\frac 14}\sqrt{\frac \pi 2}\sqrt {X^2+Y^2}, \end{eqnarray} \tag{ 28 }$$

$$\begin{eqnarray} y&=& 3^{\frac 14}\sqrt{\frac 2\pi}\sqrt {X^2+Y^2} \arctan \frac YX. \end{eqnarray} \tag{ 29 }$$

Then, in order to establish the expressions for the other sextants, it is enough to calculate $\mathcal R \circ \mathcal T_0^{-1}\circ \mathcal R^T$ and $\mathcal R^T \circ \mathcal T_0^{-1}\circ \mathcal R$, since

$$\begin{eqnarray*}\hspace{4pt}(\mathcal R \circ \mathcal T_0\circ \mathcal R^T)^{-1}=&\mathcal R \circ \mathcal T_0^{-1}\circ \mathcal R^T,\\ (\mathcal R_2 \circ \mathcal T_0\circ \mathcal R_2^T)^{-1}=&\mathcal R_2 \circ \mathcal T_0^{-1}\circ \mathcal R_2^T=-\mathcal R^T \circ \mathcal T_0^{-1}\circ \mathcal R,\\ (\mathcal R_3 \circ \mathcal T_0\circ \mathcal R_3^T)^{-1}=&\hspace{-2pt}-\mathcal T_0^{-1},\\ (\mathcal R_4 \circ \mathcal T_0\circ \mathcal R_4^T)^{-1}=&\hspace{-2pt}-\mathcal R \circ \mathcal T_0^{-1}\circ \mathcal R^T,\\ (\mathcal R_5 \circ \mathcal T_0\circ \mathcal R_5^T)^{-1}=&\mathcal R^T \circ \mathcal T_0^{-1}\circ \mathcal R. \end{eqnarray*}$$

Finally, after some calculations, we find that $\mathcal T^{-1}_H$ maps a point (X, Y)_H ≠ (0, 0) into the point (x, y), as follows:

• For (X, Y) ∈ I₀ ∪ I₃
  $$\begin{equation} \fl (x,y)=3^{-\frac 14}\sqrt {X^2+Y^2}\,{\rm sign}(X)\left( \sqrt{\frac \pi 2},\,\sqrt{\frac 6\pi} \arctan \frac YX \right); \end{equation} \tag{ 30a }$$
• For (X, Y) ∈ I₁ ∪ I₄,
  $$\begin{eqnarray} &&\fl (x,y)=\frac{3^\frac 14\sqrt {X^2+Y^2}}{\sqrt {2\pi}}\,\mbox{sign}(X)\left( \sqrt3\Big(\frac \pi 6-\arctan \frac{Y-\sqrt 3 X}{X+\sqrt 3 Y}\,\Big),\right.\nonumber\\ &&\left.\frac \pi 2+\arctan \frac{Y-\sqrt 3 X}{X+\sqrt 3 Y} \right); \end{eqnarray} \tag{ 30b }$$
• For (X, Y) ∈ I₂ ∪ I₅,
  $$\begin{eqnarray} &&\fl (x,y)=\frac{3^\frac 14\sqrt {X^2+Y^2}}{\sqrt {2\pi}}\,\mbox{sign}(X)\left( \sqrt3\Big(\frac \pi 6+\arctan \frac{Y+\sqrt 3 X}{X-\sqrt 3 Y}\,\Big),\right.\nonumber\\ &&\left.-\frac \pi 2+\arctan \frac{Y+\sqrt 3 X}{X-\sqrt 3Y} \right). \end{eqnarray} \tag{ 30c }$$

2.5. The inverse equilateral triangle map

Defining $\theta= {y\pi}/{3\sqrt 3 x}$ and $\eta=3^{\frac 34}/\sqrt{\pi}$, we obtain for the inverse mapping in sector J₀ of the equilateral triangle:

$$\begin{equation} x=\frac{\sqrt{X^2+Y^2}}\eta \tqs {\rm and}\tqs y=\frac{3\sqrt 3}{\pi}\,x \arctan \frac YX. \end{equation} \tag{ 31 }$$

Therefore, the inverse mapping for the first sector maps the point (X, Y) ∈ J₀ into the point (x, y) ∈ J₀ given by

$$\begin{eqnarray} x&=& 3^{-\frac 34}\sqrt{\pi}\sqrt {X^2+Y^2}, \end{eqnarray} \tag{ 32a }$$

$$\begin{eqnarray} y&=& \frac{3^{\frac 34}}{\sqrt{\pi}}\sqrt {X^2+Y^2} \arctan \frac YX. \end{eqnarray} \tag{ 32b }$$

Applying the rotations from the previous section for even values of k we find that $\mathcal T^{-1}_{T}$ maps a point (X, Y)_T ≠ (0, 0) into the point (x, y) as follows:

• For (X, Y) ∈ J₀
  $$\begin{equation} \fl (x,y)=3^{-\frac 34}\sqrt {X^2+Y^2}\,\left( \sqrt{\pi},\,3\sqrt{\frac 3\pi} \arctan \frac YX \right); \end{equation} \tag{ 33a }$$
• For (X, Y) ∈ J₁
  $$\begin{eqnarray} &&\fl (x,y)=-\frac{3^{\frac 14}}{2\sqrt{\pi}}\sqrt {X^2+Y^2}\,\left( \frac{\pi}{3}+3\arctan\frac{\sqrt{3}X+Y}{X-\sqrt{3}Y},\right.\nonumber\\ &&\left.-\frac{\pi}{\sqrt{3}}+\sqrt{3}\arctan\frac{\sqrt{3}X+Y}{X-\sqrt{3}Y}\right); \end{eqnarray} \tag{ 33b }$$
• For (X, Y) ∈ J₂
  $$\begin{eqnarray} &&\fl (x,y)=-\frac{3^{\frac 14}}{2\sqrt{\pi}}\sqrt {X^2+Y^2}\,\left( \frac{\pi}{3}+3\arctan\frac{\sqrt{3}X-Y}{X+\sqrt{3}Y},\right.\nonumber\\ &&\left.\frac{\pi}{\sqrt{3}}-\sqrt{3}\arctan\frac{\sqrt{3}X-Y}{X+\sqrt{3}Y}\right); \end{eqnarray} \tag{ 33c }$$

2.6. Mapping onto the sphere

Let $\mathbb S^2$ be the unit sphere centered at O and let $\mathbb S^2_\pm$ denote its northern and southern hemispheres, respectively. Let $\mathcal L_N:\mathcal C_{\sqrt 2}\to \mathbb S^2_+$ be the inverse Lambert equal-area projection with respect to the North Pole, defined by

$$\begin{equation} \fl (\widetilde{x},\widetilde{y},\widetilde{z}) = \left(\frac{X}{2}\sqrt{4-(X^2+Y^2)} ,\frac{Y}{2}\sqrt{4-(X^2+Y^2)} , 1-\frac{1}{2}(X^2+Y^2)\right) \in \mathbb S^2_+ , \end{equation} \tag{ 34 }$$

for all $(X,Y)\in \mathcal C_{\sqrt 2}$. A complete description of all known spherical projections from a sphere or parts of a sphere to the plane, used in cartography, is realized in [9]. Note that the radius of the circle $\mathcal C_R$ that is projected onto the (unit-radius) hemisphere must be equal to $R=\sqrt{2}$. For the southern hemisphere $\widetilde z$ changes into $-\widetilde z$.

The Lambert mapping $\mathcal L_N^{-1}:\mathbb S^2_+ \to \mathcal C_{\sqrt 2}$ is given by

$$\begin{equation} (X,Y) = \sqrt{\frac{2}{1+\widetilde{z}}} \left(\widetilde{x},\widetilde{y} \right)\in\mathcal{C}_{\sqrt{2}}, \end{equation} \tag{ 35 }$$

and for $\mathbb S^2_{-}$, in (35) $\widetilde z$ changes into $-\widetilde z$.

The regular hexagon H_a that maps onto the circle with radius $R=\sqrt{2}$ has the edge $a=2\sqrt \pi 3^{-3/4}$, and therefore $\alpha=3^{1/4}\sqrt{\pi}$. For the equilateral triangle we have $b=2\sqrt{2\pi}3^{-1/4}$ and $\beta=\sqrt{2\pi}3^{3/4}$. The mapping from the hexagonal grid to the northern hemisphere is then carried out by combining (25a), (25b) and (25c) with (34), and the inverse mapping can be obtained from the combination of (35) with (30a), (30b) and (30c). For the equilateral triangle T_b one combines (26a), (26b) and (26c) with (34) for the forward projection onto the northern hemisphere, and (35), (33a), (33b) and (33c) for the inverse mapping.

Figure 4 shows two grids on the sphere, the images of the planar grids in figure 3.

Zoom In Zoom Out Reset image size

Consider the regular hexagon H_a with $a=2\sqrt \pi 3^{-3/4}$. By dividing each triangle in the sextants I_k, k = 1, ..., 6, into n small triangles of equal area, we obtain a uniform grid on H_a. Mapping this grid on the two hemispheres of $\mathbb S^2$ we obtain a spherical grid having 12n² cells of equal area. The set of their vertices will be denoted by Ω_n. Similarly, for the equilateral triangle T_b, with $b=2\sqrt{2\pi}3^{-1/4}$, we construct the set Γ_n of points of $\mathbb S^2$. One has |Ω_n| = 6n² + 2 and |Γ_n| = 3n² + 2.

In the following, in order to measure the uniformity of the configurations of points Γ_n and Ω_n, we investigate their energies. For a set ω_m = {x₁, ..., x_m} of distinct points on $\mathbb S^2$ and a fixed $s\in \mathbb R$, the Riesz s-energy of the configuration ω_m is defined as

$$\begin{equation*} E_{s}(\omega_{m})=\left\{\begin{array}{@{}ll@{}} \displaystyle\sum_{1\leq i<j\leq m} \ln {\|x_{i}-x_{j}\|}^{-1},& s=0,\\ \displaystyle\sum_{1\leq i<j\leq m}\|x_{i}-x_{j}\|^{-s}, & s\neq 0, \end{array} \right . \end{equation*}$$

where ||·|| denotes the Euclidean norm.

For each point x_i ∈ ω_m, the point s-energies E_s(x_i) are defined as

$$\begin{equation*} E_{s}(x_{i})=\left\{\begin{array}{@{}ll@{}} \displaystyle \sum_{j\neq i} \ln \|x_{i}-x_{j}\|^{-1},& s=0,\\\displaystyle \sum_{j\neq i}\|x_{i}-x_{j}\|^{-s}, & s\neq 0. \end{array} \right. \end{equation*}$$

In figures 5 and 6 we have plotted the point energies for Ω₆ and Ω₁₀, and in figures 7 and 8 the point energies for Γ₈ and Γ₁₄.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In table 1 we give the values of E_max = max{E_s(x), x ∈ Ω_n} and E_min = min{E_s(x), x ∈ Ω_n}, for n = 1, ..., 15 and s = −1, 0, 1. In table 2 we give the values of E_max = max{E_s(x), x ∈ Γ_n} and E_min = min{E_s(x), x ∈ Γ_n}, for n = 1, ..., 15 and s = −1, 0, 1.

We have also compared the energies of our configurations with the extremal s-energies of m points, defined as

$$\begin{equation} \mathcal{E}_{s}(m)=\left\{\begin{array}{@{}ll@{}} \displaystyle\inf_{\omega_{m}\subset \mathbb S^{2}} E_{s}(\omega_{m}),& s\geq 0, \\ \displaystyle \sup_{\omega_{m}\subset \mathbb S^{2}} E_{s}(\omega_{m}), & s<0, \end{array} \right. \end{equation} \tag{ 36 }$$

where inf and sup are taken over all m-points subsets of $\mathbb S^2$. The determination of s-extremal configurations is a difficult problem of constraint optimization and has applications in physics, chemistry and computer science. For s = 1, the minimization of (36) is the classical Thompson problem of electrons restricted to the sphere and interacting through the Coulomb potential. This problem is relevant in molecular modeling and in the study of certain virus structures. Up to now, one could calculate extremal energies for m ⩽ 200 points and s = −1, 0, 1 and these energies are listed in [10]. We have compared the energies of our configurations Ω_n and Γ_n with the optimal ones, i.e. we have calculated the energies E_s(Ω_n) for n = 1, ..., 5 and E_s(Γ_n) for n = 1, ..., 8 (see tables 3 and 4, respectively).

Table 3. Riesz s-energies E_s(Ω_n) and extremal s-energies $\mathcal E_s(M)$ (M = 6n² + 2), for s = −1, 0, 1.

M	s = −1	s = 0	s = 1
8	E = 41.3629	E = −10.2273	E = 19.9494
	$\mathcal{E}=41.4731$	$\mathcal{E}=-10.4280$	$\mathcal{E}=19.6753$
26	E = 448.3299	E = −86.0668	E = 267.5817
	$\mathcal{E}=448.5856$	$\mathcal{E}=-87.0094$	$\mathcal{E}=265.1333$
56	E = 2087.3101	E = −358.7790	E = 1343.8734
	$\mathcal{E}=2087.6412$	$\mathcal{E}=-360.5459$	$\mathcal{E}=1337.0949$
98	E = 6398.2823	E = −1039.4824	E = 4281.1864
	$\mathcal{E}=6398.6801$	$\mathcal{E}=-1042.2847$	$\mathcal{E}=4266.8225$
152	E = 15 397.2488	E = −2421.9798	E = 10 543.8732
	$\mathcal{E}=15\,397.7081$	$\mathcal{E}=-2426.0176$	$\mathcal{E}=10\,517.8676$

Table 4. Riesz s-energies E_s(Ω_n) and extremal s-energies $\mathcal E_s(M)$ (M = 3n² + 2), for s = −1, 0, 1.

M	s = −1	s = 0	s = 1
5	E = 15.6814	E = −4.4205	E = 6.4747
	$\mathcal{E}=15.6814$	$\mathcal{E}=-4.4205$	$\mathcal{E}=6.4747$
14	E = 128.8638	E = −27.8053	E = 70.3686
	$\mathcal{E}=129.1204$	$\mathcal{E}=-28.4078$	$\mathcal{E}=69.3064$
29	E = 558.1092	E = −105.0794	E = 337.6106
	$\mathcal{E}=558.4722$	$\mathcal{E}=-106.2546$	$\mathcal{E}=334.6344$
50	E = 1663.3661	E = −289.7002	E = 1061.2942
	$\mathcal{E}=1663.8078$	$\mathcal{E}=-291.5286$	$\mathcal{E}=1055.1823$
77	E = 3948.6258	E = −655.5598	E = 2602.5717
	$\mathcal{E}=3949.1302$	$\mathcal{E}=-658.1178$	$\mathcal{E}=2591.8502$
110	E = 8061.8862	E = −1297.2073	E = 5430.5968
	$\mathcal{E}=8062.4440$	$\mathcal{E}=-1300.5711$	$\mathcal{E}=5413.5493$
149	E = 14 795.1468	E = −2329.9250	E = 10 122.5225
	$\mathcal{E}=14\,795.7571$	$\mathcal{E}=-2334.2218$	$\mathcal{E}=10\,096.8599$
194	E = 25 084.4073	E = −3889.7685	E = 17 363.5015
	$\mathcal{E}=25\,085.0685$	$\mathcal{E}=-3895.1169$	$\mathcal{E}=17\,326.6161$

The original motivation for the derivation of the mappings in the preceding section was the fact that simulations of EBSD patterns for crystal systems with hexagonal or trigonal Bravais lattices would benefit from the availability of an area-preserving mapping based on a 2D uniform hexagonal sampling grid. In such simulations [11], one computes a 'master EBSD pattern', which is a representation of the back-scattered electron (BSE) yield as a function of the direction in which the electron leaves a crystal. Typically, such a pattern is represented on a spherical surface. From the master pattern one extracts individual EBSD patterns, ideally through a simple interpolation process. The reason for desiring a uniform hexagonal planar sampling grid for hexagonal and trigonal materials is the fact that a square grid cannot properly represent the six-fold and three-fold rotational symmetries of these crystal systems. To reduce the number of simulations that need to be carried out to obtain the master EBSD pattern, one makes use of the point symmetry group of the crystal system, so that only a sub-region of the total orientation space needs to be sampled. Symmetry operators then copy the computed BSE yield to the correct equivalent locations.

We have implemented the hexagonal mapping $\mathcal{T}_H$ in combination with a dynamical electron scattering simulation. The basic approach is as follows: a back-scattered electron is generated in the near-surface region of a sample and travels in a direction described by the wave vector k (with length |k| equal to the inverse of the electron wave length λ). The probability that such an electron will leave the sample depends on both the crystal structure of the material and on the depth z inside the sample at which the back-scattering event occurs. It is common practice to represent this probability by an integral of the electron wave probability, evaluated at the scattering sites (the atoms in their lattice positions), over a depth range from the surface (z = 0) to a maximum depth z = z₀:

$$\begin{equation} \mathcal{P}({\bit k}) = \sum_{i}\frac{Z_i^2 DW_i}{z_0}\int_0^{z_0}\vert\Psi({\bit r}_i)\vert^2\,\rmd z, \end{equation} \tag{ 37 }$$

where r_i are atom positions in the unit cell, Z_i are the atomic numbers and the factors DW_i are Debye–Waller factors that represent the diffuseness of the atom cores caused by thermal vibrations. Ψ(r) is the electron wave function, which is typically represented by a Bloch wave expansion (a wave that has the periodicity of the underlying crystal lattice):

$$\begin{equation} \Psi({\bit r})=\sum_j\alpha^{(j)}\sum_{{\bit g}}C^{(j)}_{{\bit g}} {\rm e}^{2\pi\rmi ({\bit k}_{0}+(\gamma^{(j)} +\rmi q^{(j)}){\bit n}+{\bit g})\cdot{\bit r}}. \end{equation} \tag{ 38 }$$

The parameters α^(j) in (38) are the Bloch wave excitation amplitudes, $C_{{\bit g}}^{(j)}$ are the Bloch wave coefficients, the eigenvalues are represented as λ^(j) ≡ γ^(j) + iq^(j), and the vectors g are the reciprocal lattice vectors of the crystal. k₀ is the electron wave vector corrected for refraction. The total number of terms in the first summation is equal to the number of reciprocal lattice vectors g taken into account in the simulation. The Bloch wave quantities are obtained by solving a general complex eigenvalue problem, which is essentially a matrix representation of the quantum-mechanical Schrödinger equation [12]. For simplicity in this example simulation, we ignore the energy dependence of the BSE yield; more extensive details can be found in [11].

A simulation of the master EBSD pattern for pure titanium was carried out for an electron energy of 30 keV using equation (37). Figure 9(a) shows the resulting BSE yield represented on the surface of a sphere. Since the crystal structure of Ti is centrosymmetric, only the northern hemisphere is needed, and the Lambert map of this hemisphere is shown in figure 9(b); the hexagonal six-fold symmetry axis is located at the center of the Lambert projection and corresponds to the North–South axis of the sphere in (a). Applying the inverse hexagonal mapping $\mathcal{T}^{-1}_H$ then results in the hexagonal map shown in figure 9(c). While the map is represented in this figure using a standard square grid of pixels, the actual data points are located at the nodes of a uniform hexagonal sampling grid similar to that shown in figure 2(a). The total number of grid points along the horizontal direction is 1001. Note that the projection has been rotated by 90° with respect to the orientation of the hexagon in figure 1 to conform with the standard representation of crystallographic reference frames.

Zoom In Zoom Out Reset image size

In this contribution, we have derived area-preserving projections from uniform hexagonal domains and from uniform triangular domains onto circular domains. In combination with the Lambert equal-area projection from the sphere $\mathbb S^2$ to the disc, these new mappings provide direct bi-directional conversions between uniform planar grids with three-fold or six-fold rotational symmetry and a corresponding uniform grid on the sphere. Explicit expressions are given for both direct and inverse mappings between the hexagonal and triangular domains and the disc. The quality of the resulting grids on the sphere is analyzed by means of the Riesz s-energy of configurations with different numbers of sampling points; the configurational energies are also compared to extremal s-energies.

While area-preserving mappings on the sphere have many applications in science and engineering, we have selected one particular application: the representation of the channeling-modified back-scattered electron yield as a function of direction for crystal systems with hexagonal or trigonal Bravais lattices. An example simulation for hexagonal titanium at a microscope accelerating voltage of 30 kV shows how the projections can be used to represent data with hexagonal symmetry.

MDG would like to acknowledge the Office of Naval Research, contract # N00014-12-1-0075, for financial support. Daniela Roşca was supported by the Sectoral Operational Programme Human Resources Development 2007–2013 of the Romanian Ministry of Labor, Family and Social Protection through the Financial Agreement POSDRU/89/1.5/S/62557.