1 /**
2 * @file
3 * @brief math primitives
4 */
5
6 /*
7 Copyright (C) 1997-2001 Id Software, Inc.
8
9 This program is free software; you can redistribute it and/or
10 modify it under the terms of the GNU General Public License
11 as published by the Free Software Foundation; either version 2
12 of the License, or (at your option) any later version.
13
14 This program is distributed in the hope that it will be useful,
15 but WITHOUT ANY WARRANTY; without even the implied warranty of
16 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
17
18 See the GNU General Public License for more details.
19
20 You should have received a copy of the GNU General Public License
21 along with this program; if not, write to the Free Software
22 Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
23
24 */
25
26 #include "mathlib.h"
27 #include "../common/common.h"
28 #include <algorithm>
29
30 #ifndef logf
31 #define logf(x) ((float)log((double)(x)))
32 #endif
33
34 const vec2_t vec2_origin = { 0, 0 };
35 const vec3_t vec3_origin = { 0, 0, 0 };
36 const vec4_t vec4_origin = { 0, 0, 0, 0 };
37 const pos3_t pos3_origin = { 0, 0, 0 };
38
39 /** @brief cos 45 degree */
40 #define RT2 0.70710678118654752440084436210485
41
42 const GridBox GridBox::EMPTY(pos3_origin, pos3_origin);
43
44 /**
45 * @note DIRECTIONS
46 * straight
47 * 0 = x+1, y
48 * 1 = x-1, y
49 * 2 = x , y+1
50 * 3 = x , y-1
51 * diagonal
52 * 4 = x+1, y+1
53 * 5 = x-1, y-1
54 * 6 = x-1, y+1
55 * 7 = x+1, y-1
56 * @note (change in x, change in y, change in z, change in height status)
57 */
58 const vec4_t dvecs[PATHFINDING_DIRECTIONS] = {
59 { 1, 0, 0, 0}, /* E */
60 {-1, 0, 0, 0}, /* W */
61 { 0, 1, 0, 0}, /* N */
62 { 0, -1, 0, 0}, /* S */
63 { 1, 1, 0, 0}, /* NE */
64 {-1, -1, 0, 0}, /* SW */
65 {-1, 1, 0, 0}, /* NW */
66 { 1, -1, 0, 0}, /* SE */
67 { 0, 0, 1, 0}, /* CLIMB UP */
68 { 0, 0, -1, 0}, /* CLIMB DOWN */
69 { 0, 0, 0, -1}, /* STAND UP */
70 { 0, 0, 0, 1}, /* STAND DOWN */
71 { 0, 0, 0, 0}, /* UNDEFINED OPPOSITE OF FALL DOWN */
72 { 0, 0, -1, 0}, /* FALL DOWN */
73 { 0, 0, 0, 0}, /* UNDEFINED */
74 { 0, 0, 0, 0}, /* UNDEFINED */
75 { 1, 0, 1, 0}, /* UP E (Fliers only)*/
76 {-1, 0, 1, 0}, /* UP W (Fliers only) */
77 { 0, 1, 1, 0}, /* UP N (Fliers only) */
78 { 0, -1, 1, 0}, /* UP S (Fliers only) */
79 { 1, 1, 1, 0}, /* UP NE (Fliers only) */
80 {-1, -1, 1, 0}, /* UP SW (Fliers only) */
81 {-1, 1, 1, 0}, /* UP NW (Fliers only) */
82 { 1, -1, 1, 0}, /* UP SE (Fliers only) */
83 { 1, 0, 0, 0}, /* E (Fliers only)*/
84 {-1, 0, 0, 0}, /* W (Fliers only) */
85 { 0, 1, 0, 0}, /* N (Fliers only) */
86 { 0, -1, 0, 0}, /* S (Fliers only) */
87 { 1, 1, 0, 0}, /* NE (Fliers only) */
88 {-1, -1, 0, 0}, /* SW (Fliers only) */
89 {-1, 1, 0, 0}, /* NW (Fliers only) */
90 { 1, -1, 0, 0}, /* SE (Fliers only) */
91 { 1, 0, -1, 0}, /* DOWN E (Fliers only) */
92 {-1, 0, -1, 0}, /* DOWN W (Fliers only) */
93 { 0, 1, -1, 0}, /* DOWN N (Fliers only) */
94 { 0, -1, -1, 0}, /* DOWN S (Fliers only) */
95 { 1, 1, -1, 0}, /* DOWN NE (Fliers only) */
96 {-1, -1, -1, 0}, /* DOWN SW (Fliers only) */
97 {-1, 1, -1, 0}, /* DOWN NW (Fliers only) */
98 { 1, -1, -1, 0}, /* DOWN SE (Fliers only) */
99 };
100
101 /* 0:E 1:W 2:N 3:S 4:NE 5:SW 6:NW 7:SE */
102 const float dvecsn[CORE_DIRECTIONS][2] = { {1, 0}, {-1, 0}, {0, 1}, {0, -1}, {RT2, RT2}, {-RT2, -RT2}, {-RT2, RT2}, {RT2, -RT2} };
103 /** @note if you change directionAngles[PATHFINDING_DIRECTIONS], you must also change function AngleToDV */
104 /* 0:E 1: W 2:N 3:S 4:NE 5:SW 6:NW 7:SE */
105 const float directionAngles[CORE_DIRECTIONS] = { 0, 180.0f, 90.0f, 270.0f, 45.0f, 225.0f, 135.0f, 315.0f };
106
107 #define DIRECTION_EAST 0
108 #define DIRECTION_WEST 1
109 #define DIRECTION_NORTH 2
110 #define DIRECTION_SOUTH 3
111 #define DIRECTION_NORTHEAST 4
112 #define DIRECTION_SOUTHWEST 5
113 #define DIRECTION_NORTHWEST 6
114 #define DIRECTION_SOUTHEAST 7
115
116 const byte dvright[CORE_DIRECTIONS] = { DIRECTION_SOUTHEAST,
117 DIRECTION_NORTHWEST, DIRECTION_NORTHEAST, DIRECTION_SOUTHWEST,
118 DIRECTION_EAST, DIRECTION_WEST, DIRECTION_NORTH, DIRECTION_SOUTH };
119 const byte dvleft[CORE_DIRECTIONS] = { DIRECTION_NORTHEAST, DIRECTION_SOUTHWEST,
120 DIRECTION_NORTHWEST, DIRECTION_SOUTHEAST, DIRECTION_NORTH,
121 DIRECTION_SOUTH, DIRECTION_WEST, DIRECTION_EAST };
122
123
124 /**
125 * @brief Returns the index of array directionAngles[DIRECTIONS] whose value is the closest to angle
126 * @note This function allows to know the closest multiple of 45 degree of angle.
127 * @param[in] angle The angle (in degrees) which is tested.
128 * @return Corresponding index of array directionAngles[DIRECTIONS].
129 */
AngleToDir(int angle)130 int AngleToDir (int angle)
131 {
132 angle += 22;
133 /* set angle between 0 <= angle < 360 */
134 angle %= 360;
135 /* next step is because the result of angle %= 360 when angle is negative depends of the compiler
136 * (it can be between -360 < angle <= 0 or 0 <= angle < 360) */
137 if (angle < 0)
138 angle += 360;
139
140 /* get an integer quotient */
141 angle /= 45;
142
143 if (angle >= 0 && angle < CORE_DIRECTIONS) {
144 static const int anglesToDV[8] = {0, 4, 2, 6, 1, 5, 3, 7};
145 return anglesToDV[angle];
146 }
147
148 /* This is the default for unknown values. */
149 Com_Printf("Error in AngleToDV: shouldn't have reached this line\n");
150 return 0;
151 }
152
153 /**
154 * @brief Round to nearest integer
155 */
Q_rint(const vec_t in)156 vec_t Q_rint (const vec_t in)
157 {
158 /* round x down to the nearest integer */
159 return floor(in + 0.5);
160 }
161
162 /**
163 * @brief Calculate distance on the geoscape.
164 * @param[in] pos1 Position at the start.
165 * @param[in] pos2 Position at the end.
166 * @return Distance from pos1 to pos2.
167 * @note distance is an angle! This is the angle (in degrees) between the 2 vectors
168 * starting at earth's center and ending at pos1 or pos2. (if you prefer distance,
169 * this is also the distance on a globe of radius 180 / pi = 57)
170 */
GetDistanceOnGlobe(const vec2_t pos1,const vec2_t pos2)171 double GetDistanceOnGlobe (const vec2_t pos1, const vec2_t pos2)
172 {
173 /* convert into rad */
174 const double latitude1 = pos1[1] * torad;
175 const double latitude2 = pos2[1] * torad;
176 const double deltaLongitude = (pos1[0] - pos2[0]) * torad;
177 double distance;
178
179 distance = cos(latitude1) * cos(latitude2) * cos(deltaLongitude) + sin(latitude1) * sin(latitude2);
180 distance = std::min(std::max(-1.0, distance), 1.0);
181 distance = acos(distance) * todeg;
182
183 assert(distance >= 0.0);
184 return distance;
185 }
186
187 /**
188 * @brief
189 */
ColorNormalize(const vec3_t in,vec3_t out)190 vec_t ColorNormalize (const vec3_t in, vec3_t out)
191 {
192 float max;
193
194 /* find the brightest component */
195 max = in[0];
196 if (in[1] > max)
197 max = in[1];
198 if (in[2] > max)
199 max = in[2];
200
201 /* avoid FPE */
202 if (EQUAL(max, 0.0)) {
203 VectorClear(out);
204 return 0;
205 }
206
207 VectorScale(in, 1.0 / max, out);
208
209 return max;
210 }
211
212 /**
213 * @brief Checks whether the given vector @c v1 is closer to @c comp as the vector @c v2
214 * @param[in] v1 Vector to check whether it's closer to @c comp as @c v2
215 * @param[in] v2 Vector to check against
216 * @param[in] comp The vector to check the distance from
217 * @return Returns true if @c v1 is closer to @c comp as @c v2
218 */
VectorNearer(const vec3_t v1,const vec3_t v2,const vec3_t comp)219 bool VectorNearer (const vec3_t v1, const vec3_t v2, const vec3_t comp)
220 {
221 vec3_t d1, d2;
222
223 VectorSubtract(comp, v1, d1);
224 VectorSubtract(comp, v2, d2);
225
226 return VectorLength(d1) < VectorLength(d2);
227 }
228
229 /**
230 * @brief Calculated the normal vector for a given vec3_t
231 * @param[in] v Vector to normalize
232 * @param[out] out The normalized vector
233 * @sa VectorNormalize
234 * @return vector length as vec_t
235 * @sa CrossProduct
236 */
VectorNormalize2(const vec3_t v,vec3_t out)237 vec_t VectorNormalize2 (const vec3_t v, vec3_t out)
238 {
239 float length;
240
241 length = DotProduct(v, v);
242 length = sqrt(length); /** @todo */
243
244 if (!EQUAL(length, 0.0)) {
245 const float ilength = 1.0 / length;
246 out[0] = v[0] * ilength;
247 out[1] = v[1] * ilength;
248 out[2] = v[2] * ilength;
249 }
250
251 return length;
252 }
253
254 /**
255 * @brief Sets vector_out (vc) to vevtor1 (va) + scale * vector2 (vb)
256 * @param[in] veca Position to start from
257 * @param[in] scale Speed of the movement
258 * @param[in] vecb Movement direction
259 * @param[out] outVector Target vector
260 */
VectorMA(const vec3_t veca,const float scale,const vec3_t vecb,vec3_t outVector)261 void VectorMA (const vec3_t veca, const float scale, const vec3_t vecb, vec3_t outVector)
262 {
263 outVector[0] = veca[0] + scale * vecb[0];
264 outVector[1] = veca[1] + scale * vecb[1];
265 outVector[2] = veca[2] + scale * vecb[2];
266 }
267
VectorClampMA(vec3_t veca,float scale,const vec3_t vecb,vec3_t vecc)268 void VectorClampMA (vec3_t veca, float scale, const vec3_t vecb, vec3_t vecc)
269 {
270 int i;
271
272 /* clamp veca to bounds */
273 for (i = 0; i < 3; i++)
274 if (veca[i] > 4094.0)
275 veca[i] = 4094.0;
276 else if (veca[i] < -4094.0)
277 veca[i] = -4094.0;
278
279 /* rescale to fit */
280 for (i = 0; i < 3; i++) {
281 const float test = veca[i] + scale * vecb[i];
282 if (test < -4095.0f) {
283 const float newScale = (-4094.0 - veca[i]) / vecb[i];
284 if (fabs(newScale) < fabs(scale))
285 scale = newScale;
286 } else if (test > 4095.0f) {
287 const float newScale = (4094.0 - veca[i]) / vecb[i];
288 if (fabs(newScale) < fabs(scale))
289 scale = newScale;
290 }
291 }
292
293 /* use rescaled scale */
294 VectorMA(veca, scale, vecb, vecc);
295 }
296
297 /**
298 * @brief Multiply 3*3 matrix by 3*3 matrix.
299 * @param[out] c The result of the multiplication matrix = @c a * @c b (not the same as @c b * @c a !)
300 * @param[in] a First matrix.
301 * @param[in] b Second matrix.
302 * @sa GLMatrixMultiply
303 */
MatrixMultiply(const vec3_t a[3],const vec3_t b[3],vec3_t c[3])304 void MatrixMultiply (const vec3_t a[3], const vec3_t b[3], vec3_t c[3])
305 {
306 c[0][0] = a[0][0] * b[0][0] + a[1][0] * b[0][1] + a[2][0] * b[0][2];
307 c[0][1] = a[0][1] * b[0][0] + a[1][1] * b[0][1] + a[2][1] * b[0][2];
308 c[0][2] = a[0][2] * b[0][0] + a[1][2] * b[0][1] + a[2][2] * b[0][2];
309
310 c[1][0] = a[0][0] * b[1][0] + a[1][0] * b[1][1] + a[2][0] * b[1][2];
311 c[1][1] = a[0][1] * b[1][0] + a[1][1] * b[1][1] + a[2][1] * b[1][2];
312 c[1][2] = a[0][2] * b[1][0] + a[1][2] * b[1][1] + a[2][2] * b[1][2];
313
314 c[2][0] = a[0][0] * b[2][0] + a[1][0] * b[2][1] + a[2][0] * b[2][2];
315 c[2][1] = a[0][1] * b[2][0] + a[1][1] * b[2][1] + a[2][1] * b[2][2];
316 c[2][2] = a[0][2] * b[2][0] + a[1][2] * b[2][1] + a[2][2] * b[2][2];
317 }
318
319 /**
320 * @brief Builds an opengl translation and rotation matrix
321 * @param origin The translation vector
322 * @param angles The angle vector that is used to calculate the rotation part of the matrix
323 * @param matrix The resulting matrix, must be of dimension 16
324 */
GLMatrixAssemble(const vec3_t origin,const vec3_t angles,float * matrix)325 void GLMatrixAssemble (const vec3_t origin, const vec3_t angles, float* matrix)
326 {
327 /* fill in edge values */
328 matrix[3] = matrix[7] = matrix[11] = 0.0;
329 matrix[15] = 1.0;
330
331 /* add rotation */
332 AngleVectors(angles, &matrix[0], &matrix[4], &matrix[8]);
333 /* flip an axis */
334 VectorInverse(&matrix[4]);
335
336 /* add translation */
337 matrix[12] = origin[0];
338 matrix[13] = origin[1];
339 matrix[14] = origin[2];
340 }
341
342 /**
343 * @brief Multiply 4*4 matrix by 4*4 matrix.
344 * @sa MatrixMultiply
345 * @param[out] c The result of the multiplication matrix = @c a * @c b (not the same as @c b * @c a as matrix
346 * multiplication is not commutative)
347 * @param[in] a First matrix.
348 * @param[in] b Second matrix.
349 */
GLMatrixMultiply(const float a[16],const float b[16],float c[16])350 void GLMatrixMultiply (const float a[16], const float b[16], float c[16])
351 {
352 int i, j;
353
354 for (j = 0; j < 4; j++) {
355 int k = j * 4;
356 for (i = 0; i < 4; i++)
357 c[i + k] = a[i] * b[k] + a[i + 4] * b[k + 1] + a[i + 8] * b[k + 2] + a[i + 12] * b[k + 3];
358 }
359 }
360
361 /**
362 * @brief Multiply 4*4 matrix by 4d vector.
363 * @param[out] out the result of the multiplication = @c m * @c in.
364 * @param[in] m 4*4 matrix
365 * @param[in] in 4d vector.
366 * @sa VectorRotate
367 */
GLVectorTransform(const float m[16],const vec4_t in,vec4_t out)368 void GLVectorTransform (const float m[16], const vec4_t in, vec4_t out)
369 {
370 int i;
371
372 for (i = 0; i < 4; i++)
373 out[i] = m[i] * in[0] + m[i + 4] * in[1] + m[i + 8] * in[2] + m[i + 12] * in[3];
374 }
375
376 /**
377 * @brief Transform position (xyz) vector by OpenGL rules
378 * @note Equivalent to calling GLVectorTransfrom with (x y z 1) vector and taking first 3 components of result
379 * @param[out] out the result of the multiplication = @c m * @c in.
380 * @param[in] m 4*4 matrix
381 * @param[in] in 3d vector.
382 * @sa GLVectorTransform
383 */
GLPositionTransform(const float m[16],const vec3_t in,vec3_t out)384 void GLPositionTransform (const float m[16], const vec3_t in, vec3_t out)
385 {
386 int i;
387
388 for (i = 0; i < 3; i++)
389 out[i] = m[i] * in[0] + m[i + 4] * in[1] + m[i + 8] * in[2] + m[i + 12];
390 }
391
392 /**
393 * @brief Rotate a vector with a rotation matrix.
394 * @param[out] vb The result of multiplication (ie vector @c va after rotation).
395 * @param[in] m The 3*3 matrix (rotation matrix in case of a rotation).
396 * @param[in] va The vector that should be multiplied (ie rotated in case of rotation).
397 * @note Basically, this is just the multiplication of @c m * @c va: this is the same than
398 * GLVectorTransform in 3D. This can be used for other applications than rotation.
399 * @sa GLVectorTransform
400 */
VectorRotate(vec3_t m[3],const vec3_t va,vec3_t vb)401 void VectorRotate (vec3_t m[3], const vec3_t va, vec3_t vb)
402 {
403 vb[0] = m[0][0] * va[0] + m[1][0] * va[1] + m[2][0] * va[2];
404 vb[1] = m[0][1] * va[0] + m[1][1] * va[1] + m[2][1] * va[2];
405 vb[2] = m[0][2] * va[0] + m[1][2] * va[1] + m[2][2] * va[2];
406 }
407
408 /**
409 * @brief Compare two vectors that may have an epsilon difference but still be
410 * the same vectors
411 * @param[in] v1 vector to compare with v2
412 * @param[in] v2 vector to compare with v1
413 * @param[in] epsilon The epsilon the vectors may differ to still be the same
414 * @return 1 if the 2 vectors are the same (less than epsilon difference).
415 * @note This is not an exact calculation (should use a sqrt). Use this function
416 * only if you want to know if both vectors are the same with a precison that is
417 * roughly epsilon (the difference should be lower than sqrt(3) * epsilon).
418 */
VectorCompareEps(const vec3_t v1,const vec3_t v2,float epsilon)419 int VectorCompareEps (const vec3_t v1, const vec3_t v2, float epsilon)
420 {
421 vec3_t d;
422
423 VectorSubtract(v1, v2, d);
424 d[0] = fabs(d[0]);
425 d[1] = fabs(d[1]);
426 d[2] = fabs(d[2]);
427
428 if (d[0] > epsilon || d[1] > epsilon || d[2] > epsilon)
429 return 0;
430
431 return 1;
432 }
433
434 /**
435 * @brief Calculate the length of a vector
436 * @param[in] v Vector to calculate the length of
437 * @sa VectorNormalize
438 * @return Vector length as vec_t
439 */
VectorLength(const vec3_t v)440 vec_t VectorLength (const vec3_t v)
441 {
442 return sqrtf(DotProduct(v, v));
443 }
444
445 /**
446 * @brief Calculate a position on @c v1 @c v2 line.
447 * @param[in] v1 First point of the line.
448 * @param[in] v2 Second point of the line.
449 * @param[in] mix Position on the line. If 0 < @c mix < 1, @c out is between @c v1 and @c v2 .
450 * if @c mix < 0, @c out is outside @c v1 and @c v2 , on @c v1 side.
451 * @param[out] out The resulting point
452 */
VectorMix(const vec3_t v1,const vec3_t v2,float mix,vec3_t out)453 void VectorMix (const vec3_t v1, const vec3_t v2, float mix, vec3_t out)
454 {
455 const float number = 1.0 - mix;
456
457 out[0] = v1[0] * number + v2[0] * mix;
458 out[1] = v1[1] * number + v2[1] * mix;
459 out[2] = v1[2] * number + v2[2] * mix;
460 }
461
462 /**
463 * @brief Inverse a vector.
464 * @param[in,out] v Vector to inverse. Output value is -@c v.
465 */
VectorInverse(vec3_t v)466 void VectorInverse (vec3_t v)
467 {
468 v[0] = -v[0];
469 v[1] = -v[1];
470 v[2] = -v[2];
471 }
472
473 /**
474 * @brief Calculates the midpoint between two vectors.
475 * @param[in] point1 vector of first point.
476 * @param[in] point2 vector of second point.
477 * @param[out] midpoint calculated midpoint vector.
478 */
VectorMidpoint(const vec3_t point1,const vec3_t point2,vec3_t midpoint)479 void VectorMidpoint (const vec3_t point1, const vec3_t point2, vec3_t midpoint)
480 {
481 VectorAdd(point1, point2, midpoint);
482 VectorScale(midpoint, 0.5f, midpoint);
483 }
484
485 /**
486 * @brief Calculates the angle (in radians) between the two given vectors.
487 * @note Both vectors must be normalized.
488 * @return the angle in radians.
489 */
VectorAngleBetween(const vec3_t vec1,const vec3_t vec2)490 float VectorAngleBetween (const vec3_t vec1, const vec3_t vec2)
491 {
492 const float dot = DotProduct(vec1, vec2);
493 const float angle = acos(dot);
494 return angle;
495 }
496
497
Q_log2(int val)498 int Q_log2 (int val)
499 {
500 int answer = 0;
501
502 while (val >>= 1)
503 answer++;
504 return answer;
505 }
506
507 /**
508 * @brief Return random values between 0 and 1
509 * @sa crand
510 * @sa gaussrand
511 */
frand(void)512 float frand (void)
513 {
514 return (rand() & 32767) * (1.0 / 32767);
515 }
516
517
518 /**
519 * @brief Return random values between -1 and 1
520 * @sa frand
521 * @sa gaussrand
522 */
crand(void)523 float crand (void)
524 {
525 return (rand() & 32767) * (2.0 / 32767) - 1;
526 }
527
528 /**
529 * @brief generate two gaussian distributed random numbers with median at 0 and stdev of 1
530 * @param[out] gauss1 First gaussian distributed random number
531 * @param[out] gauss2 Second gaussian distributed random number
532 * @sa crand
533 * @sa frand
534 */
gaussrand(float * gauss1,float * gauss2)535 void gaussrand (float* gauss1, float* gauss2)
536 {
537 float x1, x2, w, tmp;
538
539 do {
540 x1 = crand();
541 x2 = crand();
542 w = x1 * x1 + x2 * x2;
543 } while (w >= 1.0);
544
545 tmp = -2 * logf(w);
546 w = sqrt(tmp / w);
547 *gauss1 = x1 * w;
548 *gauss2 = x2 * w;
549 }
550 /** @brief Calculates the bounding box in absolute coordinates, also for rotating objects.
551 * WARNING: do not use this for angles other than 90, 180 or 270 !! */
CalculateMinsMaxs(const vec3_t angles,const vec3_t mins,const vec3_t maxs,const vec3_t origin,vec3_t absmin,vec3_t absmax)552 void CalculateMinsMaxs (const vec3_t angles, const vec3_t mins, const vec3_t maxs, const vec3_t origin, vec3_t absmin, vec3_t absmax)
553 {
554 /* expand for rotation */
555 if (VectorNotEmpty(angles)) {
556 vec3_t minVec, maxVec, tmpMinVec, tmpMaxVec;
557 vec3_t centerVec, halfVec, newCenterVec, newHalfVec;
558 vec3_t m[3];
559
560 /* Find the center of the extents. */
561 VectorCenterFromMinsMaxs(mins, maxs, centerVec);
562
563 /* Find the half height and half width of the extents. */
564 VectorSubtract(maxs, centerVec, halfVec);
565
566 /* Rotate the center about the origin. */
567 VectorCreateRotationMatrix(angles, m);
568 VectorRotate(m, centerVec, newCenterVec);
569 VectorRotate(m, halfVec, newHalfVec);
570
571 /* Set minVec and maxVec to bound around newCenterVec at halfVec size. */
572 VectorSubtract(newCenterVec, newHalfVec, tmpMinVec);
573 VectorAdd(newCenterVec, newHalfVec, tmpMaxVec);
574
575 /* rotation may have changed min and max of the box, so adjust it */
576 minVec[0] = std::min(tmpMinVec[0], tmpMaxVec[0]);
577 minVec[1] = std::min(tmpMinVec[1], tmpMaxVec[1]);
578 minVec[2] = std::min(tmpMinVec[2], tmpMaxVec[2]);
579 maxVec[0] = std::max(tmpMinVec[0], tmpMaxVec[0]);
580 maxVec[1] = std::max(tmpMinVec[1], tmpMaxVec[1]);
581 maxVec[2] = std::max(tmpMinVec[2], tmpMaxVec[2]);
582
583 /* Adjust the absolute mins/maxs */
584 VectorAdd(origin, minVec, absmin);
585 VectorAdd(origin, maxVec, absmax);
586 } else { /* normal */
587 VectorAdd(origin, mins, absmin);
588 VectorAdd(origin, maxs, absmax);
589 }
590 }
591
592 /**
593 * @param angles The angles to calulcate the rotation matrix for
594 * @param matrix The resulting rotation matrix. The @c right part of this matrix is inversed because
595 * of the coordinate system we are using internally.
596 * @see AnglesVectors
597 * @see VectorRotatePoint
598 */
VectorCreateRotationMatrix(const vec3_t angles,vec3_t matrix[3])599 void VectorCreateRotationMatrix (const vec3_t angles, vec3_t matrix[3])
600 {
601 AngleVectors(angles, matrix[0], matrix[1], matrix[2]);
602 VectorInverse(matrix[1]);
603 }
604
605 /**
606 * @param[out] point The vector to rotate around and the location of the rotated vector
607 * @param[in] matrix The input rotation matrix
608 * @see VectorCreateRotationMatrix
609 */
VectorRotatePoint(vec3_t point,vec3_t matrix[3])610 void VectorRotatePoint (vec3_t point, vec3_t matrix[3])
611 {
612 vec3_t tvec;
613
614 VectorCopy(point, tvec);
615
616 point[0] = DotProduct(matrix[0], tvec);
617 point[1] = DotProduct(matrix[1], tvec);
618 point[2] = DotProduct(matrix[2], tvec);
619 }
620
621 /**
622 * @brief Create the rotation matrix in order to rotate something.
623 * @param[in] angles Contains the three angles PITCH, YAW and ROLL (in degree) of rotation around idle
624 * frame ({0, 1, 0}, {0, 0, 1} ,{1, 0, 0}) (in this order!)
625 * @param[out] forward return the first line of the rotation matrix.
626 * @param[out] right return the second line of the rotation matrix.
627 * @param[out] up return the third line of the rotation matrix.
628 * @note This matrix is the product of the 3 matrixes R*P*Y (in this order!), where R is the rotation matrix around {1, 0, 0} only
629 * (angle of rotation is angle[2]), P is the rotation matrix around {0, 1, 0} only, and Y is the rotation matrix around {0, 0, 1}.
630 * @note Due to z convention for Quake, the z-axis is inverted. Therefore, if you want to use this function in a direct frame, don't
631 * forget to inverse the second line (@c right).
632 * Exemple : to rotate v2 into v :
633 * AngleVectors(angles, m[0], m[1], m[2]);
634 * VectorInverse(m[1]);
635 * VectorRotate(m, v2, v);
636 * @sa RotatePointAroundVector : If you need rather "one rotation around a vector you choose" instead of "3 rotations around 3 vectors you don't choose".
637 */
AngleVectors(const vec3_t angles,vec3_t forward,vec3_t right,vec3_t up)638 void AngleVectors (const vec3_t angles, vec3_t forward, vec3_t right, vec3_t up)
639 {
640 const float anglePitch = angles[PITCH] * torad;
641 const float sp = sin(anglePitch);
642 const float cp = cos(anglePitch);
643 const float angleYaw = angles[YAW] * torad;
644 const float sy = sin(angleYaw);
645 const float cy = cos(angleYaw);
646 const float angleRoll = angles[ROLL] * torad;
647 const float sr = sin(angleRoll);
648 const float cr = cos(angleRoll);
649
650 if (forward) {
651 forward[0] = cp * cy;
652 forward[1] = cp * sy;
653 forward[2] = -sp;
654 }
655 if (right) {
656 right[0] = (-1 * sr * sp * cy + -1 * cr * -sy);
657 right[1] = (-1 * sr * sp * sy + -1 * cr * cy);
658 right[2] = -1 * sr * cp;
659 }
660 if (up) {
661 up[0] = (cr * sp * cy + -sr * -sy);
662 up[1] = (cr * sp * sy + -sr * cy);
663 up[2] = cr * cp;
664 }
665 }
666
667 /**
668 * @brief Checks whether a point is visible from a given position
669 * @param[in] origin Origin to test from
670 * @param[in] dir Direction to test into
671 * @param[in] point This is the point we want to check the visibility for
672 */
FrustumVis(const vec3_t origin,int dir,const vec3_t point)673 bool FrustumVis (const vec3_t origin, int dir, const vec3_t point)
674 {
675 /* view frustum check */
676 vec3_t delta;
677 byte dv;
678
679 delta[0] = point[0] - origin[0];
680 delta[1] = point[1] - origin[1];
681 delta[2] = 0;
682 VectorNormalizeFast(delta);
683 dv = dir & (DIRECTIONS - 1);
684
685 /* test 120 frustum (cos 60 = 0.5) */
686 if ((delta[0] * dvecsn[dv][0] + delta[1] * dvecsn[dv][1]) < 0.5)
687 return false;
688
689 return true;
690 }
691
692 /**
693 * @brief Projects a point on a plane passing through the origin
694 * @param[in] point Coordinates of the point to project
695 * @param[in] normal The normal vector of the plane
696 * @param[out] dst Coordinates of the projection on the plane
697 * @pre @c Non-null pointers and a normalized normal vector.
698 */
ProjectPointOnPlane(vec3_t dst,const vec3_t point,const vec3_t normal)699 static inline void ProjectPointOnPlane (vec3_t dst, const vec3_t point, const vec3_t normal)
700 {
701 float distance; /**< closest distance from the point to the plane */
702
703 #if 0
704 vec3_t n;
705 float inv_denom;
706 /* I added a sqrt there, otherwise this function does not work for unnormalized vector (13052007 Kracken) */
707 /* old line was inv_denom = 1.0F / DotProduct(normal, normal); */
708 inv_denom = 1.0F / sqrt(DotProduct(normal, normal));
709 #endif
710
711 distance = DotProduct(normal, point);
712 #if 0
713 n[0] = normal[0] * inv_denom;
714 n[1] = normal[1] * inv_denom;
715 n[2] = normal[2] * inv_denom;
716 #endif
717
718 dst[0] = point[0] - distance * normal[0];
719 dst[1] = point[1] - distance * normal[1];
720 dst[2] = point[2] - distance * normal[2];
721 }
722
Q_rsqrtApprox(const float number)723 static inline float Q_rsqrtApprox (const float number)
724 {
725 union
726 {
727 float f;
728 int i;
729 } t;
730 float y;
731 float x2;
732 const float threehalfs = 1.5F;
733
734 x2 = number * 0.5F;
735 t.f = number;
736 /* what the fuck? */
737 t.i = 0x5f3759df - (t.i >> 1);
738 y = t.f;
739 /* 1st iteration */
740 y = y * (threehalfs - (x2 * y * y));
741 /* 2nd iteration */
742 y = y * (threehalfs - (x2 * y * y));
743 return y;
744 }
745
746 /**
747 * @brief Calculate unit vector for a given vec3_t
748 * @param[in] v Vector to normalize
749 * @sa VectorNormalize2
750 * @return vector length as vec_t
751 */
VectorNormalize(vec3_t v)752 vec_t VectorNormalize (vec3_t v)
753 {
754 const float length = sqrt(DotProduct(v, v));
755 if (length) {
756 const float ilength = 1.0 / length;
757 v[0] *= ilength;
758 v[1] *= ilength;
759 v[2] *= ilength;
760 }
761
762 return length;
763 }
764
765 /**
766 * @brief fast vector normalize routine that does not check to make sure
767 * that length != 0, nor does it return length
768 */
VectorNormalizeFast(vec3_t v)769 void VectorNormalizeFast (vec3_t v)
770 {
771 const float ilength = Q_rsqrtApprox(DotProduct(v, v));
772 v[0] *= ilength;
773 v[1] *= ilength;
774 v[2] *= ilength;
775 }
776
777 /**
778 * @brief Finds a vector perpendicular to the source vector
779 * @param[in] src The source vector
780 * @param[out] dst A vector perpendicular to @c src
781 * @note @c dst is a perpendicular vector to @c src such that it is the closest
782 * to one of the three axis: {1,0,0}, {0,1,0} and {0,0,1} (chosen in that order
783 * in case of equality)
784 * @pre non-nullptr pointers and @c src is normalized
785 * @sa ProjectPointOnPlane
786 */
PerpendicularVector(vec3_t dst,const vec3_t src)787 void PerpendicularVector (vec3_t dst, const vec3_t src)
788 {
789 int pos;
790 int i;
791 float minelem = 1.0F;
792 vec3_t tempvec;
793
794 /* find the smallest magnitude axially aligned vector */
795 for (pos = 0, i = 0; i < 3; i++) {
796 const float a = fabs(src[i]);
797 if (a < minelem) {
798 pos = i;
799 minelem = a;
800 }
801 }
802 tempvec[0] = tempvec[1] = tempvec[2] = 0.0F;
803 tempvec[pos] = 1.0F;
804
805 /* project the point onto the plane defined by src */
806 ProjectPointOnPlane(dst, tempvec, src);
807
808 /* normalize the result */
809 VectorNormalizeFast(dst);
810 }
811
812 /**
813 * @brief binary operation on vectors in a three-dimensional space
814 * @note also known as the vector product or outer product
815 * @note It differs from the dot product in that it results in a vector
816 * rather than in a scalar
817 * @note Its main use lies in the fact that the cross product of two vectors
818 * is orthogonal to both of them.
819 * @param[in] v1 directional vector
820 * @param[in] v2 directional vector
821 * @param[out] cross output
822 * @note you have the right and forward values of an axis, their cross product will
823 * @note be a properly oriented "up" direction
824 * @sa DotProduct
825 * @sa VectorNormalize2
826 */
CrossProduct(const vec3_t v1,const vec3_t v2,vec3_t cross)827 void CrossProduct (const vec3_t v1, const vec3_t v2, vec3_t cross)
828 {
829 cross[0] = v1[1] * v2[2] - v1[2] * v2[1];
830 cross[1] = v1[2] * v2[0] - v1[0] * v2[2];
831 cross[2] = v1[0] * v2[1] - v1[1] * v2[0];
832 }
833
R_ConcatRotations(float in1[3][3],float in2[3][3],float out[3][3])834 static inline void R_ConcatRotations (float in1[3][3], float in2[3][3], float out[3][3])
835 {
836 out[0][0] = in1[0][0] * in2[0][0] + in1[0][1] * in2[1][0] + in1[0][2] * in2[2][0];
837 out[0][1] = in1[0][0] * in2[0][1] + in1[0][1] * in2[1][1] + in1[0][2] * in2[2][1];
838 out[0][2] = in1[0][0] * in2[0][2] + in1[0][1] * in2[1][2] + in1[0][2] * in2[2][2];
839 out[1][0] = in1[1][0] * in2[0][0] + in1[1][1] * in2[1][0] + in1[1][2] * in2[2][0];
840 out[1][1] = in1[1][0] * in2[0][1] + in1[1][1] * in2[1][1] + in1[1][2] * in2[2][1];
841 out[1][2] = in1[1][0] * in2[0][2] + in1[1][1] * in2[1][2] + in1[1][2] * in2[2][2];
842 out[2][0] = in1[2][0] * in2[0][0] + in1[2][1] * in2[1][0] + in1[2][2] * in2[2][0];
843 out[2][1] = in1[2][0] * in2[0][1] + in1[2][1] * in2[1][1] + in1[2][2] * in2[2][1];
844 out[2][2] = in1[2][0] * in2[0][2] + in1[2][1] * in2[1][2] + in1[2][2] * in2[2][2];
845 }
846
847 /**
848 * @brief Rotate a point around a given vector
849 * @param[in] dir The vector around which to rotate
850 * @param[in] point The point to be rotated
851 * @param[in] degrees How many degrees to rotate the point by
852 * @param[out] dst The point after rotation
853 * @note Warning: @c dst must be different from @c point (otherwise the result has no meaning)
854 * @pre @c dir must be normalized
855 */
RotatePointAroundVector(vec3_t dst,const vec3_t dir,const vec3_t point,float degrees)856 void RotatePointAroundVector (vec3_t dst, const vec3_t dir, const vec3_t point, float degrees)
857 {
858 float m[3][3];
859 float im[3][3];
860 float zrot[3][3];
861 float tmpmat[3][3];
862 float rot[3][3];
863 int i;
864 vec3_t vr, vup, vf;
865
866 vf[0] = dir[0];
867 vf[1] = dir[1];
868 vf[2] = dir[2];
869
870 PerpendicularVector(vr, dir);
871 CrossProduct(vr, vf, vup);
872
873 m[0][0] = vr[0];
874 m[1][0] = vr[1];
875 m[2][0] = vr[2];
876
877 m[0][1] = vup[0];
878 m[1][1] = vup[1];
879 m[2][1] = vup[2];
880
881 m[0][2] = vf[0];
882 m[1][2] = vf[1];
883 m[2][2] = vf[2];
884
885 memcpy(im, m, sizeof(im));
886
887 im[0][1] = m[1][0];
888 im[0][2] = m[2][0];
889 im[1][0] = m[0][1];
890 im[1][2] = m[2][1];
891 im[2][0] = m[0][2];
892 im[2][1] = m[1][2];
893
894 OBJZERO(zrot);
895
896 /* now prepare the rotation matrix */
897 zrot[0][0] = cos(degrees * torad);
898 zrot[0][1] = sin(degrees * torad);
899 zrot[1][0] = -sin(degrees * torad);
900 zrot[1][1] = cos(degrees * torad);
901 zrot[2][2] = 1.0F;
902
903 R_ConcatRotations(m, zrot, tmpmat);
904 R_ConcatRotations(tmpmat, im, rot);
905
906 for (i = 0; i < 3; i++) {
907 dst[i] = DotProduct(rot[i], point);
908 }
909 }
910
911 /**
912 * @brief Print a 3D vector
913 * @param[in] v The vector to be printed
914 * @param[out] text The resulting string. Must be pre-allocated !
915 */
Print3Vector(const vec3_t v,const char * text)916 void Print3Vector (const vec3_t v, const char* text)
917 {
918 Com_Printf("%s (%f, %f, %f)\n", text, v[0], v[1], v[2]);
919 }
920
921 /**
922 * @brief Print a 2D vector
923 * @param[in] v The vector to be printed
924 * @param[out] text The resulting string. Must be pre-allocated !
925 */
Print2Vector(const vec2_t v,const char * text)926 void Print2Vector (const vec2_t v, const char* text)
927 {
928 Com_Printf("%s (%f, %f)\n", text, v[0], v[1]);
929 }
930
931 /**
932 * @brief Converts longitude and latitude to a 3D vector in Euclidean
933 * coordinates
934 * @param[in] a The longitude and latitude in a 2D vector
935 * @param[out] v The resulting normal 3D vector
936 * @sa VecToPolar
937 */
PolarToVec(const vec2_t a,vec3_t v)938 void PolarToVec (const vec2_t a, vec3_t v)
939 {
940 const float p = a[0] * torad; /* long */
941 const float t = a[1] * torad; /* lat */
942 /* v[0] = z, v[1] = x, v[2] = y - wtf? */
943 VectorSet(v, cos(p) * cos(t), sin(p) * cos(t), sin(t));
944 }
945
946 /**
947 * @brief Converts vector coordinates into polar coordinates
948 * @sa PolarToVec
949 */
VecToPolar(const vec3_t v,vec2_t a)950 void VecToPolar (const vec3_t v, vec2_t a)
951 {
952 a[0] = todeg * atan2(v[1], v[0]); /* long */
953 a[1] = 90 - todeg * acos(v[2]); /* lat */
954 }
955
956 /**
957 * @brief Converts a vector to an angle vector
958 * @param[in] value1
959 * @param[in] angles Target vector for pitch, yaw, roll
960 * @sa anglemod
961 */
VecToAngles(const vec3_t value1,vec3_t angles)962 void VecToAngles (const vec3_t value1, vec3_t angles)
963 {
964 float yaw, pitch;
965
966 /* only check the first two values for being zero */
967 if (Vector2Empty(value1)) {
968 yaw = 0.0f;
969 if (value1[2] > 0)
970 pitch = 90.0f;
971 else
972 pitch = 270.0f;
973 } else {
974 const float forward = sqrt(value1[0] * value1[0] + value1[1] * value1[1]);
975 if (!EQUAL(value1[0], 0.0))
976 yaw = atan2(value1[1], value1[0]) * todeg;
977 else if (value1[1] > 0.0f)
978 yaw = 90.0f;
979 else
980 yaw = -90.0f;
981 if (yaw < 0.0)
982 yaw += 360.0f;
983
984 pitch = atan2(value1[2], forward) * todeg;
985 if (pitch < 0.0f)
986 pitch += 360.0f;
987 }
988
989 /* up and down */
990 angles[PITCH] = -pitch;
991 /* left and right */
992 angles[YAW] = yaw;
993 /* tilt left and right */
994 angles[ROLL] = 0.0f;
995 }
996
997 /**
998 * @brief Checks whether i is power of two value
999 */
Q_IsPowerOfTwo(int i)1000 bool Q_IsPowerOfTwo (int i)
1001 {
1002 return (i > 0 && !(i & (i - 1)));
1003 }
1004
1005 /**
1006 * @brief Returns the angle resulting from turning fraction * angle
1007 * from angle1 to angle2
1008 */
LerpAngle(float a2,float a1,float frac)1009 float LerpAngle (float a2, float a1, float frac)
1010 {
1011 if (a1 - a2 > 180)
1012 a1 -= 360;
1013 if (a1 - a2 < -180)
1014 a1 += 360;
1015 return a2 + frac * (a1 - a2);
1016 }
1017
1018 /**
1019 * @brief returns angle normalized to the range [0 <= angle < 360]
1020 * @param[in] angle Angle
1021 * @sa VecToAngles
1022 */
AngleNormalize360(float angle)1023 float AngleNormalize360 (float angle)
1024 {
1025 return (360.0 / 65536) * ((int)(angle * (65536 / 360.0)) & 65535);
1026 }
1027
1028 /**
1029 * @brief returns angle normalized to the range [-180 < angle <= 180]
1030 * @param[in] angle Angle
1031 */
AngleNormalize180(float angle)1032 float AngleNormalize180 (float angle)
1033 {
1034 angle = AngleNormalize360(angle);
1035 if (angle > 180.0)
1036 angle -= 360.0;
1037 return angle;
1038 }
1039
1040 /**
1041 * @brief Calculates the center of a bounding box
1042 * @param[in] mins The lower end of the bounding box
1043 * @param[in] maxs The upper end of the bounding box
1044 * @param[out] center The target center vector calculated from @c mins and @c maxs
1045 * @note If the mins/maxs come from an AABB, use getCenter() instead
1046 */
VectorCenterFromMinsMaxs(const vec3_t mins,const vec3_t maxs,vec3_t center)1047 void VectorCenterFromMinsMaxs (const vec3_t mins, const vec3_t maxs, vec3_t center)
1048 {
1049 VectorAdd(mins, maxs, center);
1050 VectorScale(center, 0.5, center);
1051 }
1052
1053 /**
1054 * @brief Calculates a bounding box from a center and a size
1055 * @param[in] center The center vector
1056 * @param[in] size The size vector to calculate the bbox
1057 * @param[out] mins The lower end of the bounding box
1058 * @param[out] maxs The upper end of the bounding box
1059 */
VectorCalcMinsMaxs(const vec3_t center,const vec3_t size,vec3_t mins,vec3_t maxs)1060 void VectorCalcMinsMaxs (const vec3_t center, const vec3_t size, vec3_t mins, vec3_t maxs)
1061 {
1062 int i;
1063
1064 for (i = 0; i < 3; i++) {
1065 const vec_t length = abs(size[i]) / 2;
1066 mins[i] = center[i] - length;
1067 maxs[i] = center[i] + length;
1068 }
1069 }
1070
1071 /**
1072 * @brief Sets mins and maxs to their starting points before using AddPointToBounds
1073 * @sa AddPointToBounds
1074 */
ClearBounds(vec3_t mins,vec3_t maxs)1075 void ClearBounds (vec3_t mins, vec3_t maxs)
1076 {
1077 mins[0] = mins[1] = mins[2] = 99999;
1078 maxs[0] = maxs[1] = maxs[2] = -99999;
1079 }
1080
1081 /**
1082 * @brief If the point is outside the box defined by mins and maxs, expand
1083 * the box to accommodate it. Sets mins and maxs to their new values
1084 */
AddPointToBounds(const vec3_t v,vec3_t mins,vec3_t maxs)1085 void AddPointToBounds (const vec3_t v, vec3_t mins, vec3_t maxs)
1086 {
1087 int i;
1088 for (i = 0; i < 3; i++) {
1089 vec_t val = v[i];
1090 if (val < mins[i])
1091 mins[i] = val;
1092 if (val > maxs[i])
1093 maxs[i] = val;
1094 }
1095 }
1096
1097 /**
1098 * @brief Projects the normalized directional vectors on to the normal's plane.
1099 * The fourth component of the resulting tangent vector represents sidedness.
1100 */
TangentVectors(const vec3_t normal,const vec3_t sdir,const vec3_t tdir,vec4_t tangent,vec3_t binormal)1101 void TangentVectors (const vec3_t normal, const vec3_t sdir, const vec3_t tdir, vec4_t tangent, vec3_t binormal)
1102 {
1103 vec3_t s, t;
1104
1105 /* normalize the directional vectors */
1106 VectorCopy(sdir, s);
1107 VectorNormalizeFast(s);
1108
1109 VectorCopy(tdir, t);
1110 VectorNormalizeFast(t);
1111
1112 /* project the directional vector onto the plane */
1113 VectorMA(s, -DotProduct(s, normal), normal, tangent);
1114 VectorNormalizeFast(tangent);
1115
1116 /* resolve sidedness, encode as fourth tangent component */
1117 CrossProduct(normal, tangent, binormal);
1118
1119 if (DotProduct(t, binormal) < 0.0)
1120 tangent[3] = -1.0;
1121 else
1122 tangent[3] = 1.0;
1123
1124 VectorScale(binormal, tangent[3], binormal);
1125 }
1126
1127 /**
1128 * @brief Grahm-Schmidt orthogonalization
1129 * @param[out] out Orthogonalized vector
1130 * @param[in] in Reference vector
1131 */
Orthogonalize(vec3_t out,const vec3_t in)1132 void Orthogonalize (vec3_t out, const vec3_t in)
1133 {
1134 vec3_t tmp;
1135 VectorMul(DotProduct(out, in), in, tmp);
1136 VectorSubtract(out, tmp, out);
1137 VectorNormalizeFast(out);
1138 }
1139
1140 /**
1141 * @brief Transposes @c m and stores the result in @c t
1142 * @param[in] m The matrix to transpose
1143 * @param[out] t The transposed matrix
1144 */
MatrixTranspose(const vec3_t m[3],vec3_t t[3])1145 void MatrixTranspose (const vec3_t m[3], vec3_t t[3])
1146 {
1147 int i, j;
1148
1149 for (i = 0; i < 3; i++) {
1150 for(j = 0; j < 3; j++) {
1151 t[i][j] = m[j][i];
1152 }
1153 }
1154 }
1155
RayIntersectAABB(const vec3_t start,const vec3_t end,const vec3_t mins,const vec3_t maxs)1156 bool RayIntersectAABB (const vec3_t start, const vec3_t end, const vec3_t mins, const vec3_t maxs)
1157 {
1158 float t0 = 0.0f;
1159 float t1 = 1.0f;
1160 vec3_t delta;
1161 int i;
1162
1163 VectorSubtract(end, start, delta);
1164
1165 for (i = 0; i < 3; i++) {
1166 const float threshold = 1.0e-6f;
1167 float u0, u1;
1168
1169 if (fabs(delta[i]) < threshold) {
1170 if (delta[i] > 0.0f) {
1171 return !(end[i] < mins[i] || start[i] > maxs[i]);
1172 } else {
1173 return !(start[i] < mins[i] || end[i] > maxs[i]);
1174 }
1175 }
1176
1177 u0 = (mins[i] - start[i]) / delta[i];
1178 u1 = (maxs[i] - start[i]) / delta[i];
1179
1180 if (u0 > u1) {
1181 const float temp = u0;
1182 u0 = u1;
1183 u1 = temp;
1184 }
1185
1186 if (u1 < t0 || u0 > t1) {
1187 return false;
1188 }
1189
1190 t0 = std::max(u0, t0);
1191 t1 = std::min(u1, t1);
1192
1193 if (t1 < t0) {
1194 return false;
1195 }
1196 }
1197
1198 return true;
1199 }
1200