1 /*
2 * This program is free software; you can redistribute it and/or
3 * modify it under the terms of the GNU General Public License
4 * as published by the Free Software Foundation; either version 2
5 * of the License, or (at your option) any later version.
6 *
7 * This program is distributed in the hope that it will be useful,
8 * but WITHOUT ANY WARRANTY; without even the implied warranty of
9 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
10 * GNU General Public License for more details.
11 *
12 * You should have received a copy of the GNU General Public License
13 * along with this program; if not, write to the Free Software Foundation,
14 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
15 *
16 * The Original Code is Copyright (C) 2001-2002 by NaN Holding BV.
17 * All rights reserved.
18 *
19 * The Original Code is: some of this file.
20 *
21 * */
22
23 /** \file
24 * \ingroup bli
25 */
26
27 #ifndef __MATH_GEOM_INLINE_C__
28 #define __MATH_GEOM_INLINE_C__
29
30 #include "BLI_math.h"
31 #include "BLI_math_vector.h"
32
33 #include <string.h>
34
35 /* A few small defines. Keep'em local! */
36 #define SMALL_NUMBER 1.e-8f
37
38 /********************************** Polygons *********************************/
39
cross_tri_v2(const float v1[2],const float v2[2],const float v3[2])40 MINLINE float cross_tri_v2(const float v1[2], const float v2[2], const float v3[2])
41 {
42 return (v1[0] - v2[0]) * (v2[1] - v3[1]) + (v1[1] - v2[1]) * (v3[0] - v2[0]);
43 }
44
area_tri_signed_v2(const float v1[2],const float v2[2],const float v3[2])45 MINLINE float area_tri_signed_v2(const float v1[2], const float v2[2], const float v3[2])
46 {
47 return 0.5f * ((v1[0] - v2[0]) * (v2[1] - v3[1]) + (v1[1] - v2[1]) * (v3[0] - v2[0]));
48 }
49
area_tri_v2(const float v1[2],const float v2[2],const float v3[2])50 MINLINE float area_tri_v2(const float v1[2], const float v2[2], const float v3[2])
51 {
52 return fabsf(area_tri_signed_v2(v1, v2, v3));
53 }
54
area_squared_tri_v2(const float v1[2],const float v2[2],const float v3[2])55 MINLINE float area_squared_tri_v2(const float v1[2], const float v2[2], const float v3[2])
56 {
57 float area = area_tri_signed_v2(v1, v2, v3);
58 return area * area;
59 }
60
61 /****************************** Spherical Harmonics **************************/
62
zero_sh(float r[9])63 MINLINE void zero_sh(float r[9])
64 {
65 memset(r, 0, sizeof(float[9]));
66 }
67
copy_sh_sh(float r[9],const float a[9])68 MINLINE void copy_sh_sh(float r[9], const float a[9])
69 {
70 memcpy(r, a, sizeof(float[9]));
71 }
72
mul_sh_fl(float r[9],const float f)73 MINLINE void mul_sh_fl(float r[9], const float f)
74 {
75 int i;
76
77 for (i = 0; i < 9; i++) {
78 r[i] *= f;
79 }
80 }
81
add_sh_shsh(float r[9],const float a[9],const float b[9])82 MINLINE void add_sh_shsh(float r[9], const float a[9], const float b[9])
83 {
84 int i;
85
86 for (i = 0; i < 9; i++) {
87 r[i] = a[i] + b[i];
88 }
89 }
90
dot_shsh(const float a[9],const float b[9])91 MINLINE float dot_shsh(const float a[9], const float b[9])
92 {
93 float r = 0.0f;
94 int i;
95
96 for (i = 0; i < 9; i++) {
97 r += a[i] * b[i];
98 }
99
100 return r;
101 }
102
diffuse_shv3(float sh[9],const float v[3])103 MINLINE float diffuse_shv3(float sh[9], const float v[3])
104 {
105 /* See formula (13) in:
106 * "An Efficient Representation for Irradiance Environment Maps" */
107 static const float c1 = 0.429043f, c2 = 0.511664f, c3 = 0.743125f;
108 static const float c4 = 0.886227f, c5 = 0.247708f;
109 float x, y, z, sum;
110
111 x = v[0];
112 y = v[1];
113 z = v[2];
114
115 sum = c1 * sh[8] * (x * x - y * y);
116 sum += c3 * sh[6] * z * z;
117 sum += c4 * sh[0];
118 sum += -c5 * sh[6];
119 sum += 2.0f * c1 * (sh[4] * x * y + sh[7] * x * z + sh[5] * y * z);
120 sum += 2.0f * c2 * (sh[3] * x + sh[1] * y + sh[2] * z);
121
122 return sum;
123 }
124
vec_fac_to_sh(float r[9],const float v[3],const float f)125 MINLINE void vec_fac_to_sh(float r[9], const float v[3], const float f)
126 {
127 /* See formula (3) in:
128 * "An Efficient Representation for Irradiance Environment Maps" */
129 float sh[9], x, y, z;
130
131 x = v[0];
132 y = v[1];
133 z = v[2];
134
135 sh[0] = 0.282095f;
136
137 sh[1] = 0.488603f * y;
138 sh[2] = 0.488603f * z;
139 sh[3] = 0.488603f * x;
140
141 sh[4] = 1.092548f * x * y;
142 sh[5] = 1.092548f * y * z;
143 sh[6] = 0.315392f * (3.0f * z * z - 1.0f);
144 sh[7] = 1.092548f * x * z;
145 sh[8] = 0.546274f * (x * x - y * y);
146
147 mul_sh_fl(sh, f);
148 copy_sh_sh(r, sh);
149 }
150
eval_shv3(float sh[9],const float v[3])151 MINLINE float eval_shv3(float sh[9], const float v[3])
152 {
153 float tmp[9];
154
155 vec_fac_to_sh(tmp, v, 1.0f);
156 return dot_shsh(tmp, sh);
157 }
158
madd_sh_shfl(float r[9],const float sh[9],const float f)159 MINLINE void madd_sh_shfl(float r[9], const float sh[9], const float f)
160 {
161 float tmp[9];
162
163 copy_sh_sh(tmp, sh);
164 mul_sh_fl(tmp, f);
165 add_sh_shsh(r, r, tmp);
166 }
167
168 /* get the 2 dominant axis values, 0==X, 1==Y, 2==Z */
axis_dominant_v3(int * r_axis_a,int * r_axis_b,const float axis[3])169 MINLINE void axis_dominant_v3(int *r_axis_a, int *r_axis_b, const float axis[3])
170 {
171 const float xn = fabsf(axis[0]);
172 const float yn = fabsf(axis[1]);
173 const float zn = fabsf(axis[2]);
174
175 if (zn >= xn && zn >= yn) {
176 *r_axis_a = 0;
177 *r_axis_b = 1;
178 }
179 else if (yn >= xn && yn >= zn) {
180 *r_axis_a = 0;
181 *r_axis_b = 2;
182 }
183 else {
184 *r_axis_a = 1;
185 *r_axis_b = 2;
186 }
187 }
188
189 /* same as axis_dominant_v3 but return the max value */
axis_dominant_v3_max(int * r_axis_a,int * r_axis_b,const float axis[3])190 MINLINE float axis_dominant_v3_max(int *r_axis_a, int *r_axis_b, const float axis[3])
191 {
192 const float xn = fabsf(axis[0]);
193 const float yn = fabsf(axis[1]);
194 const float zn = fabsf(axis[2]);
195
196 if (zn >= xn && zn >= yn) {
197 *r_axis_a = 0;
198 *r_axis_b = 1;
199 return zn;
200 }
201 else if (yn >= xn && yn >= zn) {
202 *r_axis_a = 0;
203 *r_axis_b = 2;
204 return yn;
205 }
206 else {
207 *r_axis_a = 1;
208 *r_axis_b = 2;
209 return xn;
210 }
211 }
212
213 /* get the single dominant axis value, 0==X, 1==Y, 2==Z */
axis_dominant_v3_single(const float vec[3])214 MINLINE int axis_dominant_v3_single(const float vec[3])
215 {
216 const float x = fabsf(vec[0]);
217 const float y = fabsf(vec[1]);
218 const float z = fabsf(vec[2]);
219 return ((x > y) ? ((x > z) ? 0 : 2) : ((y > z) ? 1 : 2));
220 }
221
222 /* the dominant axis of an orthogonal vector */
axis_dominant_v3_ortho_single(const float vec[3])223 MINLINE int axis_dominant_v3_ortho_single(const float vec[3])
224 {
225 const float x = fabsf(vec[0]);
226 const float y = fabsf(vec[1]);
227 const float z = fabsf(vec[2]);
228 return ((x < y) ? ((x < z) ? 0 : 2) : ((y < z) ? 1 : 2));
229 }
230
max_axis_v3(const float vec[3])231 MINLINE int max_axis_v3(const float vec[3])
232 {
233 const float x = vec[0];
234 const float y = vec[1];
235 const float z = vec[2];
236 return ((x > y) ? ((x > z) ? 0 : 2) : ((y > z) ? 1 : 2));
237 }
238
min_axis_v3(const float vec[3])239 MINLINE int min_axis_v3(const float vec[3])
240 {
241 const float x = vec[0];
242 const float y = vec[1];
243 const float z = vec[2];
244 return ((x < y) ? ((x < z) ? 0 : 2) : ((y < z) ? 1 : 2));
245 }
246
247 /**
248 * Simple method to find how many tri's we need when we already know the corner+poly count.
249 *
250 * \param poly_count: The number of ngon's/tris (1-2 sided faces will give incorrect results)
251 * \param corner_count: also known as loops in BMesh/DNA
252 */
poly_to_tri_count(const int poly_count,const int corner_count)253 MINLINE int poly_to_tri_count(const int poly_count, const int corner_count)
254 {
255 BLI_assert(!poly_count || corner_count > poly_count * 2);
256 return corner_count - (poly_count * 2);
257 }
258
plane_point_side_v3(const float plane[4],const float co[3])259 MINLINE float plane_point_side_v3(const float plane[4], const float co[3])
260 {
261 return dot_v3v3(co, plane) + plane[3];
262 }
263
264 /* useful to calculate an even width shell, by taking the angle between 2 planes.
265 * The return value is a scale on the offset.
266 * no angle between planes is 1.0, as the angle between the 2 planes approaches 180d
267 * the distance gets very high, 180d would be inf, but this case isn't valid */
shell_angle_to_dist(const float angle)268 MINLINE float shell_angle_to_dist(const float angle)
269 {
270 return (UNLIKELY(angle < SMALL_NUMBER)) ? 1.0f : fabsf(1.0f / cosf(angle));
271 }
272 /**
273 * equivalent to ``shell_angle_to_dist(angle_normalized_v3v3(a, b))``
274 */
shell_v3v3_normalized_to_dist(const float a[3],const float b[3])275 MINLINE float shell_v3v3_normalized_to_dist(const float a[3], const float b[3])
276 {
277 const float angle_cos = fabsf(dot_v3v3(a, b));
278 BLI_ASSERT_UNIT_V3(a);
279 BLI_ASSERT_UNIT_V3(b);
280 return (UNLIKELY(angle_cos < SMALL_NUMBER)) ? 1.0f : (1.0f / angle_cos);
281 }
282 /**
283 * equivalent to ``shell_angle_to_dist(angle_normalized_v2v2(a, b))``
284 */
shell_v2v2_normalized_to_dist(const float a[2],const float b[2])285 MINLINE float shell_v2v2_normalized_to_dist(const float a[2], const float b[2])
286 {
287 const float angle_cos = fabsf(dot_v2v2(a, b));
288 BLI_ASSERT_UNIT_V2(a);
289 BLI_ASSERT_UNIT_V2(b);
290 return (UNLIKELY(angle_cos < SMALL_NUMBER)) ? 1.0f : (1.0f / angle_cos);
291 }
292
293 /**
294 * equivalent to ``shell_angle_to_dist(angle_normalized_v3v3(a, b) / 2)``
295 */
shell_v3v3_mid_normalized_to_dist(const float a[3],const float b[3])296 MINLINE float shell_v3v3_mid_normalized_to_dist(const float a[3], const float b[3])
297 {
298 float angle_cos;
299 float ab[3];
300 BLI_ASSERT_UNIT_V3(a);
301 BLI_ASSERT_UNIT_V3(b);
302 add_v3_v3v3(ab, a, b);
303 angle_cos = (normalize_v3(ab) != 0.0f) ? fabsf(dot_v3v3(a, ab)) : 0.0f;
304 return (UNLIKELY(angle_cos < SMALL_NUMBER)) ? 1.0f : (1.0f / angle_cos);
305 }
306
307 /**
308 * equivalent to ``shell_angle_to_dist(angle_normalized_v2v2(a, b) / 2)``
309 */
shell_v2v2_mid_normalized_to_dist(const float a[2],const float b[2])310 MINLINE float shell_v2v2_mid_normalized_to_dist(const float a[2], const float b[2])
311 {
312 float angle_cos;
313 float ab[2];
314 BLI_ASSERT_UNIT_V2(a);
315 BLI_ASSERT_UNIT_V2(b);
316 add_v2_v2v2(ab, a, b);
317 angle_cos = (normalize_v2(ab) != 0.0f) ? fabsf(dot_v2v2(a, ab)) : 0.0f;
318 return (UNLIKELY(angle_cos < SMALL_NUMBER)) ? 1.0f : (1.0f / angle_cos);
319 }
320
321 #undef SMALL_NUMBER
322
323 #endif /* __MATH_GEOM_INLINE_C__ */
324