Line data Source code
1 : /**
2 : * @file test_Geometry.cpp
3 : *
4 : * @brief Test for Geometry class
5 : *
6 : * @author Tobias Lasser - initial code
7 : */
8 :
9 : #include "doctest/doctest.h"
10 : #include "Geometry.h"
11 : #include "VolumeDescriptor.h"
12 :
13 : using namespace elsa;
14 : using namespace doctest;
15 :
16 : TEST_SUITE_BEGIN("core");
17 :
18 : TEST_CASE("Geometry: Testing 2D geometries")
19 20 : {
20 20 : using namespace geometry;
21 20 : GIVEN("some 2D setup")
22 20 : {
23 20 : IndexVector_t volCoeff(2);
24 20 : volCoeff << 5, 5;
25 20 : VolumeDescriptor ddVol(volCoeff);
26 :
27 20 : IndexVector_t detCoeff(2);
28 20 : detCoeff << 5, 1;
29 20 : VolumeDescriptor ddDet(detCoeff);
30 :
31 20 : real_t s2c = 10;
32 20 : real_t c2d = 4;
33 :
34 20 : WHEN("testing geometry without rotation/offsets")
35 20 : {
36 5 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d}, Radian{0},
37 5 : VolumeData2D{Size2D{volCoeff}}, SinogramData2D{Size2D{detCoeff}});
38 :
39 5 : THEN("a copy is the same")
40 5 : {
41 1 : Geometry gCopy(g);
42 1 : REQUIRE_EQ(gCopy, g);
43 1 : }
44 :
45 5 : THEN("then P and Pinv are inverse")
46 5 : {
47 1 : RealMatrix_t id(2, 2);
48 1 : id.setIdentity();
49 :
50 1 : auto P = g.getProjectionMatrix();
51 1 : auto Pinv = g.getInverseProjectionMatrix();
52 1 : RealMatrix_t result = (P * Pinv) - id;
53 1 : REQUIRE_EQ(result.sum(), Approx(0));
54 1 : }
55 :
56 5 : THEN("then the rotation matrix is correct")
57 5 : {
58 1 : RealMatrix_t id(2, 2);
59 1 : id.setIdentity();
60 1 : REQUIRE_EQ((g.getRotationMatrix() - id).sum(), Approx(0));
61 1 : }
62 :
63 5 : THEN("the camera center is correct")
64 5 : {
65 1 : auto c = g.getCameraCenter();
66 1 : auto o = ddVol.getLocationOfOrigin();
67 1 : REQUIRE_EQ((c[0] - o[0]), Approx(0));
68 1 : REQUIRE_EQ((c[1] - o[1] + s2c), Approx(0));
69 1 : }
70 :
71 5 : THEN("the source detector distance is correct")
72 5 : {
73 1 : REQUIRE_EQ(g.getSourceDetectorDistance() - s2c - c2d, Approx(0));
74 1 : }
75 5 : }
76 :
77 20 : WHEN("testing geometry without rotation but with principal point offset")
78 20 : {
79 5 : real_t px = 2;
80 :
81 5 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d}, Radian{0},
82 5 : VolumeData2D{Size2D{volCoeff}}, SinogramData2D{Size2D{detCoeff}},
83 5 : PrincipalPointOffset{px});
84 :
85 5 : THEN("a copy is the same")
86 5 : {
87 1 : Geometry gCopy(g);
88 1 : REQUIRE_EQ(gCopy, g);
89 1 : }
90 :
91 5 : THEN("then P and Pinv are inverse")
92 5 : {
93 1 : RealMatrix_t id(2, 2);
94 1 : id.setIdentity();
95 :
96 1 : auto P = g.getProjectionMatrix();
97 1 : auto Pinv = g.getInverseProjectionMatrix();
98 1 : RealMatrix_t result = (P * Pinv) - id;
99 1 : REQUIRE_EQ(result.sum(), Approx(0));
100 1 : }
101 :
102 5 : THEN("then the rotation matrix is correct")
103 5 : {
104 1 : RealMatrix_t id(2, 2);
105 1 : id.setIdentity();
106 1 : REQUIRE_EQ((g.getRotationMatrix() - id).sum(), Approx(0));
107 1 : }
108 :
109 5 : THEN("the camera center is correct")
110 5 : {
111 1 : auto c = g.getCameraCenter();
112 1 : auto o = ddVol.getLocationOfOrigin();
113 1 : REQUIRE_EQ((c[0] - o[0]), Approx(0));
114 1 : REQUIRE_EQ((c[1] - o[1] + s2c), Approx(0));
115 1 : }
116 :
117 5 : THEN("the source detector distance is correct")
118 5 : {
119 1 : REQUIRE_EQ(g.getSourceDetectorDistance() - s2c - c2d, Approx(0));
120 1 : }
121 5 : }
122 :
123 20 : WHEN("testing geometry with 90 degree rotation and no offsets")
124 20 : {
125 5 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d}, Degree{90},
126 5 : VolumeData2D{Size2D{volCoeff}}, SinogramData2D{Size2D{detCoeff}});
127 :
128 5 : THEN("a copy is the same")
129 5 : {
130 1 : Geometry gCopy(g);
131 1 : REQUIRE_EQ(gCopy, g);
132 1 : }
133 :
134 5 : THEN("then P and Pinv are inverse")
135 5 : {
136 1 : RealMatrix_t id(2, 2);
137 1 : id.setIdentity();
138 :
139 1 : auto P = g.getProjectionMatrix();
140 1 : auto Pinv = g.getInverseProjectionMatrix();
141 1 : RealMatrix_t result = (P * Pinv) - id;
142 1 : REQUIRE_EQ(result.sum(), Approx(0.0).epsilon(0.01));
143 1 : }
144 :
145 5 : THEN("then the rotation matrix is correct")
146 5 : {
147 1 : RealMatrix_t rot(2, 2);
148 1 : rot << 0, -1, 1, 0;
149 1 : REQUIRE_EQ((g.getRotationMatrix() - rot).sum(), Approx(0).epsilon(0.01));
150 1 : }
151 :
152 5 : THEN("the camera center is correct")
153 5 : {
154 1 : auto c = g.getCameraCenter();
155 1 : auto o = ddVol.getLocationOfOrigin();
156 :
157 1 : REQUIRE_EQ((c[0] - o[0] + s2c), Approx(0));
158 1 : REQUIRE_EQ((c[1] - o[1]), Approx(0).epsilon(0.01));
159 1 : }
160 :
161 5 : THEN("the source detector distance is correct")
162 5 : {
163 1 : REQUIRE_EQ(g.getSourceDetectorDistance() - s2c - c2d, Approx(0));
164 1 : }
165 5 : }
166 :
167 20 : WHEN("testing geometry with 45 degree rotation and offset center of rotation")
168 20 : {
169 5 : real_t angle = pi_t / 4; // 45 degrees
170 5 : real_t cx = -1;
171 5 : real_t cy = 2;
172 :
173 5 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d},
174 5 : Radian{angle}, VolumeData2D{Size2D{volCoeff}},
175 5 : SinogramData2D{Size2D{detCoeff}}, PrincipalPointOffset{0},
176 5 : RotationOffset2D{cx, cy});
177 :
178 5 : THEN("a copy is the same")
179 5 : {
180 1 : Geometry gCopy(g);
181 1 : REQUIRE_EQ(gCopy, g);
182 1 : }
183 :
184 5 : THEN("then P and Pinv are inverse")
185 5 : {
186 1 : RealMatrix_t id(2, 2);
187 1 : id.setIdentity();
188 :
189 1 : auto P = g.getProjectionMatrix();
190 1 : auto Pinv = g.getInverseProjectionMatrix();
191 1 : RealMatrix_t result = (P * Pinv) - id;
192 1 : REQUIRE_EQ(result.sum(), Approx(0.0).epsilon(0.01));
193 1 : }
194 :
195 5 : THEN("then the rotation matrix is correct")
196 5 : {
197 1 : RealMatrix_t rot(2, 2);
198 1 : rot << std::cos(angle), -std::sin(angle), std::sin(angle), std::cos(angle);
199 1 : REQUIRE_EQ((g.getRotationMatrix() - rot).sum(), Approx(0));
200 1 : }
201 :
202 5 : THEN("the camera center is correct")
203 5 : {
204 1 : auto c = g.getCameraCenter();
205 1 : auto o = ddVol.getLocationOfOrigin();
206 :
207 1 : real_t oldX = 0;
208 1 : real_t oldY = -s2c;
209 1 : real_t newX = std::cos(angle) * oldX + std::sin(angle) * oldY + o[0] + cx;
210 1 : real_t newY = -std::sin(angle) * oldX + std::cos(angle) * oldY + o[1] + cy;
211 :
212 1 : REQUIRE_EQ((c[0] - newX), Approx(0).epsilon(0.01));
213 1 : REQUIRE_EQ((c[1] - newY), Approx(0).epsilon(0.01));
214 1 : }
215 :
216 5 : THEN("the source detector distance is correct")
217 5 : {
218 1 : REQUIRE_EQ(g.getSourceDetectorDistance() - s2c - c2d, Approx(0));
219 1 : }
220 5 : }
221 20 : }
222 20 : }
223 :
224 : TEST_CASE("Geometry: Testing 3D geometries")
225 30 : {
226 30 : using namespace geometry;
227 30 : GIVEN("some 3D setup")
228 30 : {
229 30 : IndexVector_t volCoeff(3);
230 30 : volCoeff << 5, 5, 5;
231 30 : VolumeDescriptor ddVol(volCoeff);
232 :
233 30 : IndexVector_t detCoeff(3);
234 30 : detCoeff << 5, 5, 1;
235 30 : VolumeDescriptor ddDet(detCoeff);
236 :
237 30 : real_t s2c = 10;
238 30 : real_t c2d = 4;
239 :
240 30 : WHEN("testing geometry without rotation/offsets")
241 30 : {
242 6 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d},
243 6 : VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
244 6 : RotationAngles3D{Radian{0}, Radian{0}, Radian{0}});
245 :
246 6 : THEN("a copy is the same")
247 6 : {
248 1 : Geometry gCopy(g);
249 1 : REQUIRE_EQ(gCopy, g);
250 1 : }
251 :
252 6 : THEN("a copy constructed with the camera and rotation/translation matrices is the same")
253 6 : {
254 1 : Geometry gNew{VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
255 1 : g.getRotationMatrix(), g.getTranslationVector(),
256 1 : g.getIntrinsicMatrix()};
257 :
258 1 : REQUIRE_EQ(gNew, g);
259 1 : }
260 :
261 6 : THEN("then P and Pinv are inverse")
262 6 : {
263 1 : RealMatrix_t id(3, 3);
264 1 : id.setIdentity();
265 :
266 1 : auto P = g.getProjectionMatrix();
267 1 : auto Pinv = g.getInverseProjectionMatrix();
268 1 : RealMatrix_t result = (P * Pinv) - id;
269 1 : REQUIRE_EQ(result.sum(), Approx(0));
270 1 : }
271 :
272 6 : THEN("then the rotation matrix is correct")
273 6 : {
274 1 : RealMatrix_t id(3, 3);
275 1 : id.setIdentity();
276 1 : REQUIRE_EQ((g.getRotationMatrix() - id).sum(), Approx(0));
277 1 : }
278 :
279 6 : THEN("the camera center is correct")
280 6 : {
281 1 : auto c = g.getCameraCenter();
282 1 : auto o = ddVol.getLocationOfOrigin();
283 1 : REQUIRE_EQ((c[0] - o[0]), Approx(0));
284 1 : REQUIRE_EQ((c[1] - o[1]), Approx(0));
285 1 : REQUIRE_EQ((c[2] - o[2] + s2c), Approx(0));
286 1 : }
287 :
288 6 : THEN("the source detector distance is correct")
289 6 : {
290 1 : REQUIRE_EQ(g.getSourceDetectorDistance() - s2c - c2d, Approx(0));
291 1 : }
292 6 : }
293 :
294 30 : WHEN("testing geometry without rotation but with principal point offsets")
295 30 : {
296 6 : real_t px = -1;
297 6 : real_t py = 3;
298 :
299 6 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d},
300 6 : VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
301 6 : RotationAngles3D{Radian{0}, Radian{0}, Radian{0}},
302 6 : PrincipalPointOffset2D{px, py});
303 :
304 6 : THEN("a copy is the same")
305 6 : {
306 1 : Geometry gCopy(g);
307 1 : REQUIRE_EQ(gCopy, g);
308 1 : }
309 :
310 6 : THEN("a copy constructed with the camera and rotation/translation matrices is the same")
311 6 : {
312 1 : Geometry gNew{VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
313 1 : g.getRotationMatrix(), g.getTranslationVector(),
314 1 : g.getIntrinsicMatrix()};
315 :
316 1 : REQUIRE_EQ(gNew, g);
317 1 : }
318 :
319 6 : THEN("then P and Pinv are inverse")
320 6 : {
321 1 : RealMatrix_t id(3, 3);
322 1 : id.setIdentity();
323 :
324 1 : auto P = g.getProjectionMatrix();
325 1 : auto Pinv = g.getInverseProjectionMatrix();
326 1 : RealMatrix_t result = (P * Pinv) - id;
327 1 : REQUIRE_EQ(result.sum(), Approx(0));
328 1 : }
329 :
330 6 : THEN("then the rotation matrix is correct")
331 6 : {
332 1 : RealMatrix_t id(3, 3);
333 1 : id.setIdentity();
334 1 : REQUIRE_EQ((g.getRotationMatrix() - id).sum(), Approx(0));
335 1 : }
336 :
337 6 : THEN("the camera center is correct")
338 6 : {
339 1 : auto c = g.getCameraCenter();
340 1 : auto o = ddVol.getLocationOfOrigin();
341 1 : REQUIRE_EQ((c[0] - o[0]), Approx(0));
342 1 : REQUIRE_EQ((c[1] - o[1]), Approx(0));
343 1 : REQUIRE_EQ((c[2] - o[2] + s2c), Approx(0));
344 1 : }
345 :
346 6 : THEN("the source detector distance is correct")
347 6 : {
348 1 : REQUIRE_EQ(g.getSourceDetectorDistance() - s2c - c2d, Approx(0));
349 1 : }
350 6 : }
351 :
352 30 : WHEN("testing geometry with 90 degree rotation and no offsets")
353 30 : {
354 6 : real_t angle = pi_t / 2;
355 :
356 6 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d},
357 6 : VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
358 6 : RotationAngles3D{Radian{angle}, Radian{0}, Radian{0}});
359 :
360 6 : THEN("a copy is the same")
361 6 : {
362 1 : Geometry gCopy(g);
363 1 : REQUIRE_EQ(gCopy, g);
364 1 : }
365 :
366 6 : THEN("a copy constructed with the camera and rotation/translation matrices is the same")
367 6 : {
368 1 : Geometry gNew{VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
369 1 : g.getRotationMatrix(), g.getTranslationVector(),
370 1 : g.getIntrinsicMatrix()};
371 :
372 1 : REQUIRE_EQ(gNew, g);
373 1 : }
374 :
375 6 : THEN("then P and Pinv are inverse")
376 6 : {
377 1 : RealMatrix_t id(3, 3);
378 1 : id.setIdentity();
379 :
380 1 : auto P = g.getProjectionMatrix();
381 1 : auto Pinv = g.getInverseProjectionMatrix();
382 1 : RealMatrix_t result = (P * Pinv) - id;
383 1 : REQUIRE_EQ(result.sum(), Approx(0).epsilon(0.01));
384 1 : }
385 :
386 6 : THEN("then the rotation matrix is correct")
387 6 : {
388 1 : RealMatrix_t rot(3, 3);
389 1 : rot << std::cos(angle), 0, std::sin(angle), 0, 1, 0, -std::sin(angle), 0,
390 1 : std::cos(angle);
391 1 : REQUIRE_EQ((g.getRotationMatrix() - rot).sum(), Approx(0));
392 1 : }
393 :
394 6 : THEN("the camera center is correct")
395 6 : {
396 1 : auto c = g.getCameraCenter();
397 1 : auto o = ddVol.getLocationOfOrigin();
398 :
399 1 : REQUIRE_EQ((c[0] - o[0] - s2c), Approx(0));
400 1 : REQUIRE_EQ((c[1] - o[1]), Approx(0).epsilon(0.01));
401 1 : REQUIRE_EQ((c[2] - o[2]), Approx(0).epsilon(0.01));
402 1 : }
403 :
404 6 : THEN("the source detector distance is correct")
405 6 : {
406 1 : REQUIRE_EQ(g.getSourceDetectorDistance() - s2c - c2d, Approx(0));
407 1 : }
408 6 : }
409 :
410 30 : WHEN("testing geometry with 45/22.5 degree rotation and offset center of rotation")
411 30 : {
412 6 : real_t angle1 = pi_t / 4;
413 6 : real_t angle2 = pi_t / 2;
414 6 : RealVector_t offset(3);
415 6 : offset << 1, -2, -1;
416 :
417 6 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d},
418 6 : VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
419 6 : RotationAngles3D{Radian{angle1}, Radian{angle2}, Radian{0}},
420 6 : PrincipalPointOffset2D{0, 0}, RotationOffset3D{offset});
421 :
422 6 : THEN("a copy is the same")
423 6 : {
424 1 : Geometry gCopy(g);
425 1 : REQUIRE_EQ(gCopy, g);
426 1 : }
427 :
428 6 : THEN("a copy constructed with the camera and rotation/translation matrices is the same")
429 6 : {
430 1 : Geometry gNew{VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
431 1 : g.getRotationMatrix(), g.getTranslationVector(),
432 1 : g.getIntrinsicMatrix()};
433 :
434 1 : REQUIRE_EQ(gNew, g);
435 1 : }
436 :
437 6 : THEN("then P and Pinv are inverse")
438 6 : {
439 1 : RealMatrix_t id(3, 3);
440 1 : id.setIdentity();
441 :
442 1 : auto P = g.getProjectionMatrix();
443 1 : auto Pinv = g.getInverseProjectionMatrix();
444 1 : RealMatrix_t result = (P * Pinv) - id;
445 1 : REQUIRE_EQ(result.sum(), Approx(0).epsilon(0.01));
446 1 : }
447 :
448 6 : THEN("then the rotation matrix is correct")
449 6 : {
450 1 : RealMatrix_t rot1(3, 3);
451 1 : rot1 << std::cos(angle1), 0, std::sin(angle1), 0, 1, 0, -std::sin(angle1), 0,
452 1 : std::cos(angle1);
453 1 : RealMatrix_t rot2(3, 3);
454 1 : rot2 << std::cos(angle2), -std::sin(angle2), 0, std::sin(angle2), std::cos(angle2),
455 1 : 0, 0, 0, 1;
456 1 : REQUIRE_EQ((g.getRotationMatrix() - rot1 * rot2).sum(), Approx(0));
457 1 : }
458 :
459 6 : THEN("the camera center is correct")
460 6 : {
461 1 : auto c = g.getCameraCenter();
462 1 : auto o = ddVol.getLocationOfOrigin();
463 :
464 1 : RealVector_t src(3);
465 1 : src[0] = 0;
466 1 : src[1] = 0;
467 1 : src[2] = -s2c;
468 1 : RealVector_t rotSrc = g.getRotationMatrix().transpose() * src + o + offset;
469 :
470 1 : REQUIRE_EQ((rotSrc - c).sum(), Approx(0).epsilon(0.01));
471 1 : }
472 :
473 6 : THEN("the source detector distance is correct")
474 6 : {
475 1 : REQUIRE_EQ(g.getSourceDetectorDistance() - s2c - c2d, Approx(0));
476 1 : }
477 6 : }
478 :
479 30 : WHEN("testing geometry with 45/22.5/12.25 degree rotation as a rotation matrix")
480 30 : {
481 6 : real_t angle1 = pi_t / 4;
482 6 : real_t angle2 = pi_t / 2;
483 6 : real_t angle3 = pi_t / 8;
484 6 : RealMatrix_t rot1(3, 3);
485 6 : rot1 << std::cos(angle1), 0, std::sin(angle1), 0, 1, 0, -std::sin(angle1), 0,
486 6 : std::cos(angle1);
487 6 : RealMatrix_t rot2(3, 3);
488 6 : rot2 << std::cos(angle2), -std::sin(angle2), 0, std::sin(angle2), std::cos(angle2), 0,
489 6 : 0, 0, 1;
490 6 : RealMatrix_t rot3(3, 3);
491 6 : rot3 << std::cos(angle3), 0, std::sin(angle3), 0, 1, 0, -std::sin(angle3), 0,
492 6 : std::cos(angle3);
493 6 : RealMatrix_t R = rot1 * rot2 * rot3;
494 :
495 6 : Geometry g(s2c, c2d, ddVol, ddDet, R);
496 :
497 6 : THEN("a copy is the same")
498 6 : {
499 1 : Geometry gCopy(g);
500 1 : REQUIRE_EQ(gCopy, g);
501 1 : }
502 :
503 6 : THEN("a copy constructed with the camera and rotation/translation matrices is the same")
504 6 : {
505 1 : Geometry gNew{VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
506 1 : g.getRotationMatrix(), g.getTranslationVector(),
507 1 : g.getIntrinsicMatrix()};
508 :
509 1 : REQUIRE_EQ(gNew, g);
510 1 : }
511 :
512 6 : THEN("then P and Pinv are inverse")
513 6 : {
514 1 : RealMatrix_t id(3, 3);
515 1 : id.setIdentity();
516 :
517 1 : auto P = g.getProjectionMatrix();
518 1 : auto Pinv = g.getInverseProjectionMatrix();
519 1 : RealMatrix_t result = (P * Pinv) - id;
520 1 : REQUIRE_EQ(result.sum(), Approx(0).epsilon(0.01));
521 1 : }
522 :
523 6 : THEN("then the rotation matrix is correct")
524 6 : {
525 1 : REQUIRE_EQ((g.getRotationMatrix() - R).sum(), Approx(0));
526 1 : }
527 :
528 6 : THEN("the camera center is correct")
529 6 : {
530 1 : auto c = g.getCameraCenter();
531 1 : auto o = ddVol.getLocationOfOrigin();
532 :
533 1 : RealVector_t src(3);
534 1 : src[0] = 0;
535 1 : src[1] = 0;
536 1 : src[2] = -s2c;
537 1 : RealVector_t rotSrc = g.getRotationMatrix().transpose() * src + o;
538 :
539 1 : REQUIRE_EQ((rotSrc - c).sum(), Approx(0).epsilon(0.01));
540 1 : }
541 :
542 6 : THEN("the source detector distance is correct")
543 6 : {
544 1 : REQUIRE_EQ(g.getSourceDetectorDistance() - s2c - c2d, Approx(0));
545 1 : }
546 6 : }
547 30 : }
548 30 : }
549 :
550 : TEST_SUITE_END();
|