Theoretica
A C++ numerical and automatic mathematical library
Loading...
Searching...
No Matches
transform.h
Go to the documentation of this file.
1
5
6#ifndef THEORETICA_TRANSFORM_H
7#define THEORETICA_TRANSFORM_H
8
9#include "./algebra.h"
10#include "../core/error.h"
11
12
13namespace theoretica {
14
15 namespace algebra {
16
17
21 template<typename Matrix>
22 inline Matrix identity(unsigned int rows = 0, unsigned int cols = 0) {
23
24 Matrix m;
25 if(rows && cols)
26 m.resize(rows, cols);
27
28 make_identity(m);
29 return m;
30 }
31
32
39 template<typename Vector, typename Matrix>
40 inline Matrix& diagonal(Matrix& res, const Vector& v) {
41
42 if(v.size() != res.cols()) {
43 TH_MATH_ERROR("algebra::mat_diagonal", v.size(), INVALID_ARGUMENT);
44 mat_error(res);
45 return res;
46 }
47
48 if(v.size() != res.rows()) {
49 TH_MATH_ERROR("algebra::mat_diagonal", v.size(), INVALID_ARGUMENT);
50 mat_error(res);
51 return res;
52 }
53
54 for (unsigned int i = 0; i < res.rows(); ++i)
55 for (unsigned int j = 0; j < res.cols(); ++j)
56 res(i, j) = (i == j) ? v[i] : 0;
57
58 return res;
59 }
60
61
70 template<typename Matrix, typename Vector>
71 inline Matrix diagonal(
72 const Vector& v, unsigned int rows = 0, unsigned int cols = 0) {
73
74 Matrix m;
75 if(rows && cols)
76 m.resize(rows, cols);
77
78 diagonal(m, v);
79 return m;
80 }
81
82
88 template<typename Matrix, typename Vector>
89 inline Matrix translation(
90 const Vector& v, unsigned int rows = 0, unsigned int cols = 0) {
91
92 Matrix m;
93 if(rows && cols)
94 m.resize(rows, cols);
95
96 if(v.size() != (m.rows() - 1)) {
97 TH_MATH_ERROR("algebra::translation", v.size(), INVALID_ARGUMENT);
98 mat_error(m);
99 return m;
100 }
101
102 make_identity(m);
103
104 // The translation matrix in projective
105 // geometry is constructed by setting
106 // the last column's values to the vector,
107 // while other elements form the identity
108 for (unsigned int i = 0; i < m.rows() - 1; ++i)
109 m(i, m.cols() - 1) = v[i];
110
111 return m;
112 }
113
114
118 template<typename Matrix>
119 inline Matrix rotation_2d(
120 real theta, unsigned int rows = 0, unsigned int cols = 0) {
121
122 Matrix m;
123 if(rows && cols)
124 m.resize(rows, cols);
125
126 if(m.rows() < 2) {
127 TH_MATH_ERROR("algebra::rotation_2d", m.rows(), INVALID_ARGUMENT);
128 mat_error(m);
129 return m;
130 }
131
132 if(m.cols() < 2) {
133 TH_MATH_ERROR("algebra::rotation_2d", m.cols(), INVALID_ARGUMENT);
134 mat_error(m);
135 return m;
136 }
137
138 const real s = sin(theta);
139 const real c = cos(theta);
140
141 if(m.rows() > 2 || m.cols() > 2)
142 make_identity(m);
143
144 m(0, 0) = c;
145 m(0, 1) = -s;
146 m(1, 0) = s;
147 m(1, 1) = c;
148
149 return m;
150 }
151
152
159 template<typename Matrix, typename Vector>
160 inline Matrix rotation_3d(
161 real theta, const Vector& axis, unsigned int rows = 0, unsigned int cols = 0) {
162
163 Matrix m;
164 if(rows && cols)
165 m.resize(rows, cols);
166
167 if(axis.size() < 3) {
168 TH_MATH_ERROR("algebra::rotation_3d", axis.size(), INVALID_ARGUMENT);
169 mat_error(m);
170 return m;
171 }
172
173 if(m.rows() < 3) {
174 TH_MATH_ERROR("algebra::rotation_3d", m.rows(), INVALID_ARGUMENT);
175 mat_error(m);
176 return m;
177 }
178
179 if(m.cols() < 3) {
180 TH_MATH_ERROR("algebra::rotation_3d", m.cols(), INVALID_ARGUMENT);
181 mat_error(m);
182 return m;
183 }
184
185 if(m.rows() > 3 || m.cols() > 3)
186 make_identity(m);
187
188 const real s = sin(theta);
189 const real c = cos(theta);
190
191 const real Rx = (real) axis[0];
192 const real Ry = (real) axis[1];
193 const real Rz = (real) axis[2];
194
195 const real cm1 = (1 - c);
196
197 m(0, 0) = c + Rx * Rx * cm1;
198 m(0, 1) = Rx * Ry * cm1 - Rz * s;
199 m(0, 2) = Rx * Rz * cm1 - Ry * s;
200
201 m(1, 0) = Ry * Rx * cm1 + Rz * s;
202 m(1, 1) = c + Ry * Ry * cm1;
203 m(1, 2) = Ry * Rz * cm1 - Rx * s;
204
205 m(2, 0) = Rz * Rx * cm1 - Ry * s;
206 m(2, 1) = Rz * Ry * cm1 + Rx * s;
207 m(2, 2) = c + Rz * Rz * cm1;
208
209 return m;
210 }
211
212
218 template<typename Matrix>
219 inline Matrix rotation_3d_xaxis(
220 real theta, unsigned int rows = 0, unsigned int cols = 0) {
221
222 Matrix m;
223 if(rows && cols)
224 m.resize(rows, cols);
225
226 if(m.rows() < 3) {
227 TH_MATH_ERROR("algebra::rotation_3d_xaxis", m.rows(), INVALID_ARGUMENT);
228 mat_error(m);
229 return m;
230 }
231
232 if(m.cols() < 3) {
233 TH_MATH_ERROR("algebra::rotation_3d_xaxis", m.cols(), INVALID_ARGUMENT);
234 mat_error(m);
235 return m;
236 }
237
238 if(m.rows() > 3 || m.cols() > 3)
239 make_identity(m);
240
241 const real s = sin(theta);
242 const real c = cos(theta);
243
244 m(0, 0) = 1;
245 m(1, 1) = c;
246 m(2, 2) = c;
247
248 m(1, 2) = -s;
249 m(2, 1) = s;
250
251 return m;
252 }
253
254
260 template<typename Matrix>
261 inline Matrix rotation_3d_yaxis(
262 real theta, unsigned int rows = 0, unsigned int cols = 0) {
263
264 Matrix m;
265 if(rows && cols)
266 m.resize(rows, cols);
267
268 if(m.rows() < 3) {
269 TH_MATH_ERROR("algebra::rotation_3d_yaxis", m.rows(), INVALID_ARGUMENT);
270 mat_error(m);
271 return m;
272 }
273
274 if(m.cols() < 3) {
275 TH_MATH_ERROR("algebra::rotation_3d_yaxis", m.cols(), INVALID_ARGUMENT);
276 mat_error(m);
277 return m;
278 }
279
280 if(m.rows() > 3 || m.cols() > 3)
281 make_identity(m);
282
283 const real s = sin(theta);
284 const real c = cos(theta);
285
286 m(0, 0) = c;
287 m(1, 1) = 1;
288 m(2, 2) = c;
289 m(3, 3) = 1;
290
291 m(0, 2) = s;
292 m(2, 0) = -s;
293
294 return m;
295 }
296
297
303 template<typename Matrix>
304 inline Matrix rotation_3d_zaxis(
305 real theta, unsigned int rows = 0, unsigned int cols = 0) {
306
307 Matrix m;
308 if(rows && cols)
309 m.resize(rows, cols);
310
311 if(m.rows() < 3) {
312 TH_MATH_ERROR("algebra::rotation_3d_zaxis", m.rows(), INVALID_ARGUMENT);
313 mat_error(m);
314 return m;
315 }
316
317 if(m.cols() < 3) {
318 TH_MATH_ERROR("algebra::rotation_3d_zaxis", m.cols(), INVALID_ARGUMENT);
319 mat_error(m);
320 return m;
321 }
322
323 if(m.rows() > 3 || m.cols() > 3)
324 make_identity(m);
325
326 const real s = sin(theta);
327 const real c = cos(theta);
328
329 m(0, 0) = c;
330 m(1, 1) = c;
331 m(2, 2) = 1;
332 m(3, 3) = 1;
333
334 m(0, 1) = -s;
335 m(1, 0) = s;
336
337 return m;
338 }
339
340
357 template<typename Matrix>
358 inline Matrix perspective(
359 real left, real right, real bottom,
360 real top, real near, real far,
361 unsigned int rows = 0, unsigned int cols = 0) {
362
363 Matrix m;
364 if(rows && cols)
365 m.resize(rows, cols);
366
367 if(m.rows() < 4) {
368 TH_MATH_ERROR("algebra::perspective", m.rows(), INVALID_ARGUMENT);
369 mat_error(m);
370 return m;
371 }
372
373 if(m.cols() < 4) {
374 TH_MATH_ERROR("algebra::perspective", m.cols(), INVALID_ARGUMENT);
375 mat_error(m);
376 return m;
377 }
378
379 mat_zeroes(m);
380
381 m(0, 0) = 2 * near / (right - left);
382 m(2, 0) = (right + left) / (right - left);
383 m(1, 1) = 2 * near / (top - bottom);
384 m(2, 1) = (top + bottom) / (top - bottom);
385 m(2, 2) = -(far + near) / (far - near);
386 m(3, 2) = -(2 * far * near) / (far - near);
387 m(2, 3) = -1;
388 m(3, 3) = 0;
389
390 return m;
391 }
392
393
396 template<typename Matrix>
397 inline Matrix perspective_fov(
398 real fov, real aspect, real near, real far,
399 unsigned int rows = 0, unsigned int cols = 0) {
400
401 Matrix m;
402 if(rows && cols)
403 m.resize(rows, cols);
404
405 if(m.rows() < 4) {
406 TH_MATH_ERROR("algebra::perspective_fov", m.rows(), INVALID_ARGUMENT);
407 mat_error(m);
408 return m;
409 }
410
411 if(m.cols() < 4) {
412 TH_MATH_ERROR("algebra::perspective_fov", m.cols(), INVALID_ARGUMENT);
413 mat_error(m);
414 return m;
415 }
416
417 const real height = near * tan(radians(fov / 2.f));
418 const real width = height * aspect;
419
420 return perspective<Matrix>(-width, width, -height, height, near, far);
421 }
422
423
425 template<typename Matrix>
426 inline Matrix ortho(
427 real left, real right, real bottom,
428 real top, real near, real far,
429 unsigned int rows = 0, unsigned int cols = 0) {
430
431 Matrix m;
432 if(rows && cols)
433 m.resize(rows, cols);
434
435 if(m.rows() < 4) {
436 TH_MATH_ERROR("algebra::ortho", m.rows(), INVALID_ARGUMENT);
437 mat_error(m);
438 return m;
439 }
440
441 if(m.cols() < 4) {
442 TH_MATH_ERROR("algebra::ortho", m.cols(), INVALID_ARGUMENT);
443 mat_error(m);
444 return m;
445 }
446
447 mat_zeroes(m);
448
449 m(0, 0) = 2 / (right - left);
450 m(3, 0) = -(right + left) / (right - left);
451 m(1, 1) = 2 / (top - bottom);
452 m(3, 1) = -(top + bottom) / (top - bottom);
453 m(2, 2) = -2 / (far - near);
454 m(3, 2) = -(far + near) / (far - near);
455
456 return m;
457 }
458
459
462 template<
463 typename Matrix, typename Vector1,
464 typename Vector2, typename Vector3>
465 inline Matrix look_at(Vector1 camera, Vector2 target, Vector3 up) {
466
467 // Construct an orthonormal basis
468
469 Vector1 x_axis, y_axis, z_axis;
470 x_axis.resize(3);
471 y_axis.resize(3);
472 z_axis.resize(3);
473
474 // z = target - camera
475 vec_diff(z_axis, target, camera);
476 make_normalized(z_axis);
477
478 // x = z X up
479 vec_copy(x_axis, cross(z_axis, up));
480 make_normalized(x_axis);
481
482 // y = x X z
483 vec_copy(y_axis, cross(x_axis, z_axis));
484
485 // Negate z_axis to have a right-handed system
486 vec_scalmul(-1.0, z_axis);
487
488 // Construct the rotation and translation matrix
489 Matrix m;
490 m.resize(4, 4);
491
492 m(0, 0) = x_axis[1];
493 m(0, 1) = x_axis[2];
494 m(0, 2) = x_axis[3];
495 m(0, 3) = dot(camera, vec_scalmul(-1.0, x_axis));
496
497 m(1, 0) = y_axis[1];
498 m(1, 1) = y_axis[2];
499 m(1, 2) = y_axis[3];
500 m(1, 3) = dot(camera, vec_scalmul(-1.0, y_axis));
501
502 m(2, 0) = z_axis[1];
503 m(2, 1) = z_axis[2];
504 m(2, 2) = z_axis[3];
505 m(2, 3) = dot(camera, vec_scalmul(-1.0, z_axis));
506
507 m(3, 0) = 0;
508 m(3, 1) = 0;
509 m(3, 2) = 0;
510 m(3, 3) = 1;
511
512 return m;
513 }
514
515
517 template<typename Matrix>
518 inline Matrix symplectic(unsigned int rows = 0, unsigned int cols = 0) {
519
520 Matrix m;
521 if(rows && cols)
522 m.resize(rows, cols);
523
524 if(rows != cols || (rows % 2 != 0)) {
525 TH_MATH_ERROR("algebra::symplectic", rows, INVALID_ARGUMENT);
526 mat_error(m);
527 return m;
528 }
529
530 const unsigned int half = m.rows() / 2;
531 mat_zeroes(m);
532
533 for (unsigned int i = 0; i < half; ++i) {
534 m(i, i + half) = 1;
535 m(i + half, i) = -1;
536 }
537
538 return m;
539 }
540
541
550 template<typename Matrix>
551 inline Matrix hilbert(unsigned int rows = 0) {
552
553 Matrix H;
554 if (rows)
555 H.resize(rows, rows);
556
557 for (unsigned int i = 0; i < H.rows(); ++i)
558 for (unsigned int j = 0; j < H.cols(); ++j)
559 H(i, j) = 1.0 / (i + j + 1);
560
561 return H;
562 }
563
564
570 template<typename Vector1, typename Vector2>
571 inline Vector1 sphere_inversion(
572 const Vector1& p, const Vector2& c = Vector2(0), real r = 1) {
573
574 Vector1 q = p - c;
575 return c + q * square(r / q.norm());
576 }
577 }
578}
579
580#endif
algebra.h
Linear algebra routines.
TH_MATH_ERROR
#define TH_MATH_ERROR(F_NAME, VALUE, EXCEPTION)
TH_MATH_ERROR is a macro which throws exceptions or modifies errno (depending on which compiling opti...
Definition error.h:219
theoretica::algebra::symplectic
Matrix symplectic(unsigned int rows=0, unsigned int cols=0)
Generate a NxN symplectic matrix where N is even.
Definition transform.h:518
theoretica::algebra::cross
Vector1 cross(const Vector1 &v1, const Vector2 &v2)
Compute the cross product between two tridimensional vectors.
Definition algebra.h:373
theoretica::algebra::look_at
Matrix look_at(Vector1 camera, Vector2 target, Vector3 up)
Return a transformation matrix that points the field of view towards a given point from the <camera> ...
Definition transform.h:465
theoretica::algebra::make_identity
Matrix & make_identity(Matrix &m)
Overwrite a matrix with the identity matrix.
Definition algebra.h:77
theoretica::algebra::rotation_2d
Matrix rotation_2d(real theta, unsigned int rows=0, unsigned int cols=0)
Returns a matrix representing a 2D rotation.
Definition transform.h:119
theoretica::algebra::make_normalized
Vector & make_normalized(Vector &v)
Normalize a given vector overwriting it.
Definition algebra.h:327
theoretica::algebra::vec_copy
Vector1 & vec_copy(Vector1 &dest, const Vector2 &src)
Copy a vector by overwriting another.
Definition algebra.h:143
theoretica::algebra::sphere_inversion
Vector1 sphere_inversion(const Vector1 &p, const Vector2 &c=Vector2(0), real r=1)
Sphere inversion of a point with respect to a sphere of radius r centered at a point c.
Definition transform.h:571
theoretica::algebra::hilbert
Matrix hilbert(unsigned int rows=0)
Construct the Hilbert matrix of arbitrary dimension.
Definition transform.h:551
theoretica::algebra::mat_error
Matrix & mat_error(Matrix &m)
Overwrite the given matrix with the error matrix with NaN values on the diagonal and zeroes everywher...
Definition algebra.h:44
theoretica::algebra::dot
auto dot(const Vector1 &v, const Vector2 &w)
Computes the dot product between two vectors, using the Hermitian form if needed.
Definition algebra.h:351
theoretica::algebra::rotation_3d_xaxis
Matrix rotation_3d_xaxis(real theta, unsigned int rows=0, unsigned int cols=0)
Returns a matrix representing a 3D rotation around the x axis.
Definition transform.h:219
theoretica::algebra::mat_zeroes
Matrix & mat_zeroes(Matrix &m)
Overwrite a matrix with all zeroes.
Definition algebra.h:93
theoretica::algebra::vec_diff
Vector2 & vec_diff(Vector1 &v1, const Vector2 &v2)
Subtract two vectors and store the result in the first vector.
Definition algebra.h:1183
theoretica::algebra::perspective_fov
Matrix perspective_fov(real fov, real aspect, real near, real far, unsigned int rows=0, unsigned int cols=0)
Returns a perspective projection matrix, using the Field of View as parameter.
Definition transform.h:397
theoretica::algebra::rotation_3d_yaxis
Matrix rotation_3d_yaxis(real theta, unsigned int rows=0, unsigned int cols=0)
Returns a matrix representing a 3D rotation around the y axis.
Definition transform.h:261
theoretica::algebra::rotation_3d
Matrix rotation_3d(real theta, const Vector &axis, unsigned int rows=0, unsigned int cols=0)
Returns a matrix representing a 3D rotation around a given axis.
Definition transform.h:160
theoretica::algebra::diagonal
Matrix & diagonal(Matrix &res, const Vector &v)
Construct a matrix with the given vector as diagonal and zeroes everywhere else, overwriting the give...
Definition transform.h:40
theoretica::algebra::ortho
Matrix ortho(real left, real right, real bottom, real top, real near, real far, unsigned int rows=0, unsigned int cols=0)
Returns an orthogonal projection matrix.
Definition transform.h:426
theoretica::algebra::identity
Matrix identity(unsigned int rows=0, unsigned int cols=0)
Returns the identity matrix.
Definition transform.h:22
theoretica::algebra::rotation_3d_zaxis
Matrix rotation_3d_zaxis(real theta, unsigned int rows=0, unsigned int cols=0)
Returns a matrix representing a 3D rotation around the z axis.
Definition transform.h:304
theoretica::algebra::perspective
Matrix perspective(real left, real right, real bottom, real top, real near, real far, unsigned int rows=0, unsigned int cols=0)
Returns a perspective projection matrix with adjustable view volume boundaries.
Definition transform.h:358
theoretica::algebra::translation
Matrix translation(const Vector &v, unsigned int rows=0, unsigned int cols=0)
Returns a translation matrix, that is, an NxN matrix which describes a translation in N-1 dimensions.
Definition transform.h:89
theoretica::algebra::vec_scalmul
Vector & vec_scalmul(Field a, Vector &v)
Multiply a vector by a scalar of any compatible type.
Definition algebra.h:634
theoretica
Main namespace of the library which contains all functions and objects.
Definition algebra.h:27
theoretica::real
double real
A real number, defined as a floating point type.
Definition constants.h:198
theoretica::radians
TH_CONSTEXPR real radians(real degrees)
Convert degrees to radians.
Definition real_analysis.h:1252
theoretica::vector_element_t
std::remove_reference_t< decltype(std::declval< Structure >()[0])> vector_element_t
Extract the type of a vector (or any indexable container) from its operator[].
Definition core_traits.h:134
theoretica::cos
dual2 cos(dual2 x)
Compute the cosine of a second order dual number.
Definition dual2_functions.h:86
theoretica::tan
dual2 tan(dual2 x)
Compute the tangent of a second order dual number.
Definition dual2_functions.h:100
theoretica::sin
dual2 sin(dual2 x)
Compute the sine of a second order dual number.
Definition dual2_functions.h:72
theoretica::square
dual2 square(dual2 x)
Return the square of a second order dual number.
Definition dual2_functions.h:23