https://mooseframework.inl.gov
Loading...
Searching...
No Matches
Euler2RGB.C
Go to the documentation of this file.
1//* This file is part of the MOOSE framework
2//* https://mooseframework.inl.gov
3//*
4//* All rights reserved, see COPYRIGHT for full restrictions
5//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6//*
7//* Licensed under LGPL 2.1, please see LICENSE for details
8//* https://www.gnu.org/licenses/lgpl-2.1.html
9
10#include "Euler2RGB.h"
11#include "MathUtils.h"
12#include "libmesh/utility.h"
13
45Point
46euler2RGB(unsigned int sd, Real phi1, Real PHI, Real phi2, unsigned int phase, unsigned int sym)
47{
48 // Define Constants
49 const Real pi = libMesh::pi;
50 const Real pi_x2 = 2.0 * pi;
51 const Real a = std::sqrt(3.0) / 2.0;
52
53 // Preallocate and zero variables
54 unsigned int index = 0;
55 unsigned int nsym = 1;
56
57 Real blue = 0.0;
58 Real chi = 0.0;
59 Real chi_min = 0.0;
60 Real chi_max = 0.0;
61 Real chi_max2 = 0.0;
62 Real eta = 0.0;
63 Real eta_min = 0.0;
64 Real eta_max = 0.0;
65 Real green = 0.0;
66 Real maxRGB = 0.0;
67 Real red = 0.0;
68
69 unsigned int ref_dir[3] = {0, 0, 0};
70 Real g[3][3] = {{0.0, 0.0, 0.0}, {0.0, 0.0, 0.0}, {0.0, 0.0, 0.0}};
71 Real hkl[3] = {0.0, 0.0, 0.0};
72 Real S[3][3] = {{0.0, 0.0, 0.0}, {0.0, 0.0, 0.0}, {0.0, 0.0, 0.0}};
73 const Real(*SymOps)[3][3];
74
75 Point RGB;
76
77 // Assign reference sample direction
78 switch (sd)
79 {
80 case 1: // 100
81 ref_dir[0] = 1;
82 ref_dir[1] = 0;
83 ref_dir[2] = 0;
84 break;
85
86 case 2: // 010
87 ref_dir[0] = 0;
88 ref_dir[1] = 1;
89 ref_dir[2] = 0;
90 break;
91
92 case 3: // 001
93 ref_dir[0] = 0;
94 ref_dir[1] = 0;
95 ref_dir[2] = 1;
96 break;
97 };
98
99 // Define symmetry operators for each of the seven crystal systems
100 const Real SymOpsCubic[24][3][3] = {
101 {{1, 0, 0}, {0, 1, 0}, {0, 0, 1}}, {{0, 0, 1}, {1, 0, 0}, {0, 1, 0}},
102 {{0, 1, 0}, {0, 0, 1}, {1, 0, 0}}, {{0, -1, 0}, {0, 0, 1}, {-1, 0, 0}},
103 {{0, -1, 0}, {0, 0, -1}, {1, 0, 0}}, {{0, 1, 0}, {0, 0, -1}, {-1, 0, 0}},
104 {{0, 0, -1}, {1, 0, 0}, {0, -1, 0}}, {{0, 0, -1}, {-1, 0, 0}, {0, 1, 0}},
105 {{0, 0, 1}, {-1, 0, 0}, {0, -1, 0}}, {{-1, 0, 0}, {0, 1, 0}, {0, 0, -1}},
106 {{-1, 0, 0}, {0, -1, 0}, {0, 0, 1}}, {{1, 0, 0}, {0, -1, 0}, {0, 0, -1}},
107 {{0, 0, -1}, {0, -1, 0}, {-1, 0, 0}}, {{0, 0, 1}, {0, -1, 0}, {1, 0, 0}},
108 {{0, 0, 1}, {0, 1, 0}, {-1, 0, 0}}, {{0, 0, -1}, {0, 1, 0}, {1, 0, 0}},
109 {{-1, 0, 0}, {0, 0, -1}, {0, -1, 0}}, {{1, 0, 0}, {0, 0, -1}, {0, 1, 0}},
110 {{1, 0, 0}, {0, 0, 1}, {0, -1, 0}}, {{-1, 0, 0}, {0, 0, 1}, {0, 1, 0}},
111 {{0, -1, 0}, {-1, 0, 0}, {0, 0, -1}}, {{0, 1, 0}, {-1, 0, 0}, {0, 0, -1}},
112 {{0, 1, 0}, {1, 0, 0}, {0, 0, -1}}, {{0, -1, 0}, {1, 0, 0}, {0, 0, 1}}};
113
114 const Real SymOpsHexagonal[12][3][3] = {{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}},
115 {{-0.5, a, 0}, {-a, -0.5, 0}, {0, 0, 1}},
116 {{-0.5, -a, 0}, {a, -0.5, 0}, {0, 0, 1}},
117 {{0.5, a, 0}, {-a, 0.5, 0}, {0, 0, 1}},
118 {{-1, 0, 0}, {0, -1, 0}, {0, 0, 1}},
119 {{0.5, -a, 0}, {a, 0.5, 0}, {0, 0, 1}},
120 {{-0.5, -a, 0}, {-a, 0.5, 0}, {0, 0, -1}},
121 {{1, 0, 0}, {0, -1, 0}, {0, 0, -1}},
122 {{-0.5, a, 0}, {a, 0.5, 0}, {0, 0, -1}},
123 {{0.5, a, 0}, {a, -0.5, 0}, {0, 0, -1}},
124 {{-1, 0, 0}, {0, 1, 0}, {0, 0, -1}},
125 {{0.5, -a, 0}, {-a, -0.5, 0}, {0, 0, -1}}};
126
127 const Real SymOpsTetragonal[8][3][3] = {{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}},
128 {{-1, 0, 0}, {0, 1, 0}, {0, 0, -1}},
129 {{1, 0, 0}, {0, -1, 0}, {0, 0, -1}},
130 {{-1, 0, 0}, {0, -1, 0}, {0, 0, 1}},
131 {{0, 1, 0}, {-1, 0, 0}, {0, 0, 1}},
132 {{0, -1, 0}, {1, 0, 0}, {0, 0, 1}},
133 {{0, 1, 0}, {1, 0, 0}, {0, 0, -1}},
134 {{0, -1, 0}, {-1, 0, 0}, {0, 0, -1}}};
135
136 const Real SymOpsTrigonal[6][3][3] = {{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}},
137 {{-0.5, a, 0}, {-a, -0.5, 0}, {0, 0, 1}},
138 {{-0.5, -a, 0}, {a, -0.5, 0}, {0, 0, 1}},
139 {{0.5, a, 0}, {a, -0.5, 0}, {0, 0, -1}},
140 {{-1, 0, 0}, {0, 1, 0}, {0, 0, 1}},
141 {{0.5, -a, 0}, {-a, -0.5, 0}, {0, 0, -1}}};
142
143 const Real SymOpsOrthorhombic[4][3][3] = {{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}},
144 {{-1, 0, 0}, {0, 1, 0}, {0, 0, -1}},
145 {{1, 0, 0}, {0, -1, 0}, {0, 0, 1}},
146 {{-1, 0, 0}, {0, -1, 0}, {0, 0, 1}}};
147
148 const Real SymOpsMonoclinic[2][3][3] = {{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}},
149 {{-1, 0, 0}, {0, 1, 0}, {0, 0, -1}}};
150
151 const Real SymOpsTriclinic[1][3][3] = {{{1, 0, 0}, {0, 1, 0}, {0, 0, 1}}};
152
153 // Assign parameters based on crystal class (sym)
154 // Load cubic parameters (class 432)
155 if (sym == 43)
156 {
157 nsym = 24;
158 SymOps = SymOpsCubic;
159 eta_min = 0 * (pi / 180);
160 eta_max = 45 * (pi / 180);
161 chi_min = 0 * (pi / 180);
162 chi_max = std::acos(std::sqrt(1.0 / (2.0 + (Utility::pow<2>(std::tan(eta_max))))));
163 }
164
165 // Load hexagonal parameters (class 622)
166 else if (sym == 62)
167 {
168 nsym = 12;
169 SymOps = SymOpsHexagonal;
170 eta_min = 0 * (pi / 180);
171 eta_max = 30 * (pi / 180);
172 chi_min = 0 * (pi / 180);
173 chi_max = pi / 2;
174 }
175
176 // Load tetragonal parameters (class 422)
177 else if (sym == 42)
178 {
179 nsym = 8;
180 SymOps = SymOpsTetragonal;
181 eta_min = 0 * (pi / 180);
182 eta_max = 45 * (pi / 180);
183 chi_min = 0 * (pi / 180);
184 chi_max = pi / 2;
185 }
186
187 // Load trigonal parameters (class 32)
188 else if (sym == 32)
189 {
190 nsym = 6;
191 SymOps = SymOpsTrigonal;
192 eta_min = 0 * (pi / 180);
193 eta_max = 60 * (pi / 180);
194 chi_min = 0 * (pi / 180);
195 chi_max = pi / 2;
196 }
197
198 // Load orthorhombic parameters (class 22)
199 else if (sym == 22)
200 {
201 nsym = 4;
202 SymOps = SymOpsOrthorhombic;
203 eta_min = 0 * (pi / 180);
204 eta_max = 90 * (pi / 180);
205 chi_min = 0 * (pi / 180);
206 chi_max = pi / 2;
207 }
208
209 // Load monoclinic parameters (class 2)
210 else if (sym == 2)
211 {
212 nsym = 2;
213 SymOps = SymOpsMonoclinic;
214 eta_min = 0 * (pi / 180);
215 eta_max = 180 * (pi / 180);
216 chi_min = 0 * (pi / 180);
217 chi_max = pi / 2;
218 }
219
220 // Load triclinic parameters (class 1)
221 else if (sym == 1)
222 {
223 nsym = 1;
224 SymOps = SymOpsTriclinic;
225 eta_min = 0 * (pi / 180);
226 eta_max = 360 * (pi / 180);
227 chi_min = 0 * (pi / 180);
228 chi_max = pi / 2;
229 }
230
231 // Accomodate non-conforming (bad) data points
232 else
233 {
234 nsym = 0;
235 }
236
237 // Start of main routine //
238 // Assign black RGB values for bad data points (nsym = 0) or voids (phase = 0)
239 if (nsym == 0 || phase == 0)
240 RGB = 0;
241
242 // Assign black RGB value for Euler angles outside of allowable range
243 else if (phi1 > pi_x2 || PHI > pi || phi2 > pi_x2)
244 RGB = 0;
245
246 // Routine for valid set of Euler angles
247 else
248 {
249 // Construct 3X3 orientation matrix
250 g[0][0] = std::cos(phi1) * std::cos(phi2) - std::sin(phi1) * std::cos(PHI) * std::sin(phi2);
251 g[0][1] = std::sin(phi1) * std::cos(phi2) + std::cos(phi1) * std::cos(PHI) * std::sin(phi2);
252 g[0][2] = std::sin(phi2) * std::sin(PHI);
253 g[1][0] = -std::cos(phi1) * std::sin(phi2) - std::sin(phi1) * std::cos(PHI) * std::cos(phi2);
254 g[1][1] = -std::sin(phi1) * std::sin(phi2) + std::cos(phi1) * std::cos(PHI) * std::cos(phi2);
255 g[1][2] = std::cos(phi2) * std::sin(PHI);
256 g[2][0] = std::sin(phi1) * std::sin(PHI);
257 g[2][1] = -std::cos(phi1) * std::sin(PHI);
258 g[2][2] = std::cos(PHI);
259
260 // Loop to sort Euler angles into standard stereographic triangle (SST)
261 index = 0;
262 while (index < nsym)
263 {
264 // Form orientation matrix
265 for (unsigned int i = 0; i < 3; ++i)
266 for (unsigned int j = 0; j < 3; ++j)
267 {
268 S[i][j] = 0.0;
269 for (unsigned int k = 0; k < 3; ++k)
270 S[i][j] += SymOps[index][i][k] * g[k][j];
271 }
272
273 // Multiple orientation matrix by reference sample direction
274 for (unsigned int i = 0; i < 3; ++i)
275 {
276 hkl[i] = 0;
277 for (unsigned int j = 0; j < 3; ++j)
278 hkl[i] += S[i][j] * ref_dir[j];
279 }
280
281 // Convert to spherical coordinates (ignore "r" variable since r=1)
282 eta = std::abs(std::atan2(hkl[1], hkl[0]));
283 chi = std::acos(std::abs(hkl[2]));
284
285 // Continue if eta and chi values are within the SST
286 if (eta >= eta_min && eta < eta_max && chi >= chi_min && chi < chi_max)
287 break;
288
289 // Increment to next symmetry operator if not in SST
290 else
291 index++;
292
293 // Check for solution
294 mooseAssert(index != nsym, "Euler2RGB failed to map the supplied Euler angle into the SST!");
295 }
296
297 // Adjust maximum chi value to ensure it falls within the SST (cubic materials only)
298 if (sym == 43)
299 chi_max2 = std::acos(std::sqrt(1.0 / (2.0 + (Utility::pow<2>(std::tan(eta))))));
300 else
301 chi_max2 = pi / 2;
302
303 // Calculate the RGB color values and make adjustments to maximize colorspace
304 red = std::abs(1.0 - (chi / chi_max2));
305 blue = std::abs((eta - eta_min) / (eta_max - eta_min));
306 green = 1.0 - blue;
307
308 blue = blue * (chi / chi_max2);
309 green = green * (chi / chi_max2);
310
311 // Check for negative RGB values before taking square root
312 mooseAssert(red >= 0 || green >= 0 || blue >= 0, "RGB component values must be positive!");
313
314 RGB(0) = std::sqrt(red);
315 RGB(1) = std::sqrt(green);
316 RGB(2) = std::sqrt(blue);
317
318 // Find maximum value of red, green, or blue
319 maxRGB = std::max({RGB(0), RGB(1), RGB(2)});
320
321 // Normalize position of SST center point
322 RGB /= maxRGB;
323 }
324
325 return RGB;
326}
euler2RGB
Point euler2RGB(unsigned int sd, Real phi1, Real PHI, Real phi2, unsigned int phase, unsigned int sym)
This function rotates a set of three Bunge Euler angles into the standard Stereographic triangle,...
Definition Euler2RGB.C:46
eta
Point eta
libMesh::pi
const Real pi
Generated on Fri Aug 21 2026 13:41:26 for https://mooseframework.inl.gov by doxygen 1.9.8