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 16 : {
20 16 : using namespace geometry;
21 16 : GIVEN("some 2D setup")
22 16 : {
23 16 : IndexVector_t volCoeff(2);
24 16 : volCoeff << 5, 5;
25 16 : VolumeDescriptor ddVol(volCoeff);
26 :
27 16 : IndexVector_t detCoeff(2);
28 16 : detCoeff << 5, 1;
29 16 : VolumeDescriptor ddDet(detCoeff);
30 :
31 16 : real_t s2c = 10;
32 16 : real_t c2d = 4;
33 :
34 16 : WHEN("testing geometry without rotation/offsets")
35 16 : {
36 4 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d}, Radian{0},
37 4 : VolumeData2D{Size2D{volCoeff}}, SinogramData2D{Size2D{detCoeff}});
38 :
39 4 : THEN("a copy is the same")
40 4 : {
41 1 : Geometry gCopy(g);
42 1 : REQUIRE_EQ(gCopy, g);
43 1 : }
44 :
45 4 : THEN("then P and Pinv are inverse")
46 4 : {
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 4 : THEN("then the rotation matrix is correct")
57 4 : {
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 4 : THEN("the camera center is correct")
64 4 : {
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 4 : }
71 :
72 16 : WHEN("testing geometry without rotation but with principal point offset")
73 16 : {
74 4 : real_t px = 2;
75 :
76 4 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d}, Radian{0},
77 4 : VolumeData2D{Size2D{volCoeff}}, SinogramData2D{Size2D{detCoeff}},
78 4 : PrincipalPointOffset{px});
79 :
80 4 : THEN("a copy is the same")
81 4 : {
82 1 : Geometry gCopy(g);
83 1 : REQUIRE_EQ(gCopy, g);
84 1 : }
85 :
86 4 : THEN("then P and Pinv are inverse")
87 4 : {
88 1 : RealMatrix_t id(2, 2);
89 1 : id.setIdentity();
90 :
91 1 : auto P = g.getProjectionMatrix();
92 1 : auto Pinv = g.getInverseProjectionMatrix();
93 1 : RealMatrix_t result = (P * Pinv) - id;
94 1 : REQUIRE_EQ(result.sum(), Approx(0));
95 1 : }
96 :
97 4 : THEN("then the rotation matrix is correct")
98 4 : {
99 1 : RealMatrix_t id(2, 2);
100 1 : id.setIdentity();
101 1 : REQUIRE_EQ((g.getRotationMatrix() - id).sum(), Approx(0));
102 1 : }
103 :
104 4 : THEN("the camera center is correct")
105 4 : {
106 1 : auto c = g.getCameraCenter();
107 1 : auto o = ddVol.getLocationOfOrigin();
108 1 : REQUIRE_EQ((c[0] - o[0]), Approx(0));
109 1 : REQUIRE_EQ((c[1] - o[1] + s2c), Approx(0));
110 1 : }
111 4 : }
112 :
113 16 : WHEN("testing geometry with 90 degree rotation and no offsets")
114 16 : {
115 4 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d}, Degree{90},
116 4 : VolumeData2D{Size2D{volCoeff}}, SinogramData2D{Size2D{detCoeff}});
117 :
118 4 : THEN("a copy is the same")
119 4 : {
120 1 : Geometry gCopy(g);
121 1 : REQUIRE_EQ(gCopy, g);
122 1 : }
123 :
124 4 : THEN("then P and Pinv are inverse")
125 4 : {
126 1 : RealMatrix_t id(2, 2);
127 1 : id.setIdentity();
128 :
129 1 : auto P = g.getProjectionMatrix();
130 1 : auto Pinv = g.getInverseProjectionMatrix();
131 1 : RealMatrix_t result = (P * Pinv) - id;
132 1 : REQUIRE_EQ(result.sum(), Approx(0.0).epsilon(0.01));
133 1 : }
134 :
135 4 : THEN("then the rotation matrix is correct")
136 4 : {
137 1 : RealMatrix_t rot(2, 2);
138 1 : rot << 0, -1, 1, 0;
139 1 : REQUIRE_EQ((g.getRotationMatrix() - rot).sum(), Approx(0).epsilon(0.01));
140 1 : }
141 :
142 4 : THEN("the camera center is correct")
143 4 : {
144 1 : auto c = g.getCameraCenter();
145 1 : auto o = ddVol.getLocationOfOrigin();
146 :
147 1 : REQUIRE_EQ((c[0] - o[0] + s2c), Approx(0));
148 1 : REQUIRE_EQ((c[1] - o[1]), Approx(0).epsilon(0.01));
149 1 : }
150 4 : }
151 :
152 16 : WHEN("testing geometry with 45 degree rotation and offset center of rotation")
153 16 : {
154 4 : real_t angle = pi_t / 4; // 45 degrees
155 4 : real_t cx = -1;
156 4 : real_t cy = 2;
157 :
158 4 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d},
159 4 : Radian{angle}, VolumeData2D{Size2D{volCoeff}},
160 4 : SinogramData2D{Size2D{detCoeff}}, PrincipalPointOffset{0},
161 4 : RotationOffset2D{cx, cy});
162 :
163 4 : THEN("a copy is the same")
164 4 : {
165 1 : Geometry gCopy(g);
166 1 : REQUIRE_EQ(gCopy, g);
167 1 : }
168 :
169 4 : THEN("then P and Pinv are inverse")
170 4 : {
171 1 : RealMatrix_t id(2, 2);
172 1 : id.setIdentity();
173 :
174 1 : auto P = g.getProjectionMatrix();
175 1 : auto Pinv = g.getInverseProjectionMatrix();
176 1 : RealMatrix_t result = (P * Pinv) - id;
177 1 : REQUIRE_EQ(result.sum(), Approx(0.0).epsilon(0.01));
178 1 : }
179 :
180 4 : THEN("then the rotation matrix is correct")
181 4 : {
182 1 : RealMatrix_t rot(2, 2);
183 1 : rot << std::cos(angle), -std::sin(angle), std::sin(angle), std::cos(angle);
184 1 : REQUIRE_EQ((g.getRotationMatrix() - rot).sum(), Approx(0));
185 1 : }
186 :
187 4 : THEN("the camera center is correct")
188 4 : {
189 1 : auto c = g.getCameraCenter();
190 1 : auto o = ddVol.getLocationOfOrigin();
191 :
192 1 : real_t oldX = 0;
193 1 : real_t oldY = -s2c;
194 1 : real_t newX = std::cos(angle) * oldX + std::sin(angle) * oldY + o[0] + cx;
195 1 : real_t newY = -std::sin(angle) * oldX + std::cos(angle) * oldY + o[1] + cy;
196 :
197 1 : REQUIRE_EQ((c[0] - newX), Approx(0).epsilon(0.01));
198 1 : REQUIRE_EQ((c[1] - newY), Approx(0).epsilon(0.01));
199 1 : }
200 4 : }
201 16 : }
202 16 : }
203 :
204 : TEST_CASE("Geometry: Testing 3D geometries")
205 20 : {
206 20 : using namespace geometry;
207 20 : GIVEN("some 3D setup")
208 20 : {
209 20 : IndexVector_t volCoeff(3);
210 20 : volCoeff << 5, 5, 5;
211 20 : VolumeDescriptor ddVol(volCoeff);
212 :
213 20 : IndexVector_t detCoeff(3);
214 20 : detCoeff << 5, 5, 1;
215 20 : VolumeDescriptor ddDet(detCoeff);
216 :
217 20 : real_t s2c = 10;
218 20 : real_t c2d = 4;
219 :
220 20 : WHEN("testing geometry without rotation/offsets")
221 20 : {
222 4 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d},
223 4 : VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
224 4 : RotationAngles3D{Radian{0}, Radian{0}, Radian{0}});
225 :
226 4 : THEN("a copy is the same")
227 4 : {
228 1 : Geometry gCopy(g);
229 1 : REQUIRE_EQ(gCopy, g);
230 1 : }
231 :
232 4 : THEN("then P and Pinv are inverse")
233 4 : {
234 1 : RealMatrix_t id(3, 3);
235 1 : id.setIdentity();
236 :
237 1 : auto P = g.getProjectionMatrix();
238 1 : auto Pinv = g.getInverseProjectionMatrix();
239 1 : RealMatrix_t result = (P * Pinv) - id;
240 1 : REQUIRE_EQ(result.sum(), Approx(0));
241 1 : }
242 :
243 4 : THEN("then the rotation matrix is correct")
244 4 : {
245 1 : RealMatrix_t id(3, 3);
246 1 : id.setIdentity();
247 1 : REQUIRE_EQ((g.getRotationMatrix() - id).sum(), Approx(0));
248 1 : }
249 :
250 4 : THEN("the camera center is correct")
251 4 : {
252 1 : auto c = g.getCameraCenter();
253 1 : auto o = ddVol.getLocationOfOrigin();
254 1 : REQUIRE_EQ((c[0] - o[0]), Approx(0));
255 1 : REQUIRE_EQ((c[1] - o[1]), Approx(0));
256 1 : REQUIRE_EQ((c[2] - o[2] + s2c), Approx(0));
257 1 : }
258 4 : }
259 :
260 20 : WHEN("testing geometry without rotation but with principal point offsets")
261 20 : {
262 4 : real_t px = -1;
263 4 : real_t py = 3;
264 :
265 4 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d},
266 4 : VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
267 4 : RotationAngles3D{Radian{0}, Radian{0}, Radian{0}},
268 4 : PrincipalPointOffset2D{px, py});
269 :
270 4 : THEN("a copy is the same")
271 4 : {
272 1 : Geometry gCopy(g);
273 1 : REQUIRE_EQ(gCopy, g);
274 1 : }
275 :
276 4 : THEN("then P and Pinv are inverse")
277 4 : {
278 1 : RealMatrix_t id(3, 3);
279 1 : id.setIdentity();
280 :
281 1 : auto P = g.getProjectionMatrix();
282 1 : auto Pinv = g.getInverseProjectionMatrix();
283 1 : RealMatrix_t result = (P * Pinv) - id;
284 1 : REQUIRE_EQ(result.sum(), Approx(0));
285 1 : }
286 :
287 4 : THEN("then the rotation matrix is correct")
288 4 : {
289 1 : RealMatrix_t id(3, 3);
290 1 : id.setIdentity();
291 1 : REQUIRE_EQ((g.getRotationMatrix() - id).sum(), Approx(0));
292 1 : }
293 :
294 4 : THEN("the camera center is correct")
295 4 : {
296 1 : auto c = g.getCameraCenter();
297 1 : auto o = ddVol.getLocationOfOrigin();
298 1 : REQUIRE_EQ((c[0] - o[0]), Approx(0));
299 1 : REQUIRE_EQ((c[1] - o[1]), Approx(0));
300 1 : REQUIRE_EQ((c[2] - o[2] + s2c), Approx(0));
301 1 : }
302 4 : }
303 :
304 20 : WHEN("testing geometry with 90 degree rotation and no offsets")
305 20 : {
306 4 : real_t angle = pi_t / 2;
307 :
308 4 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d},
309 4 : VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
310 4 : RotationAngles3D{Radian{angle}, Radian{0}, Radian{0}});
311 :
312 4 : THEN("a copy is the same")
313 4 : {
314 1 : Geometry gCopy(g);
315 1 : REQUIRE_EQ(gCopy, g);
316 1 : }
317 :
318 4 : THEN("then P and Pinv are inverse")
319 4 : {
320 1 : RealMatrix_t id(3, 3);
321 1 : id.setIdentity();
322 :
323 1 : auto P = g.getProjectionMatrix();
324 1 : auto Pinv = g.getInverseProjectionMatrix();
325 1 : RealMatrix_t result = (P * Pinv) - id;
326 1 : REQUIRE_EQ(result.sum(), Approx(0).epsilon(0.01));
327 1 : }
328 :
329 4 : THEN("then the rotation matrix is correct")
330 4 : {
331 1 : RealMatrix_t rot(3, 3);
332 1 : rot << std::cos(angle), 0, std::sin(angle), 0, 1, 0, -std::sin(angle), 0,
333 1 : std::cos(angle);
334 1 : REQUIRE_EQ((g.getRotationMatrix() - rot).sum(), Approx(0));
335 1 : }
336 :
337 4 : THEN("the camera center is correct")
338 4 : {
339 1 : auto c = g.getCameraCenter();
340 1 : auto o = ddVol.getLocationOfOrigin();
341 :
342 1 : REQUIRE_EQ((c[0] - o[0] - s2c), Approx(0));
343 1 : REQUIRE_EQ((c[1] - o[1]), Approx(0).epsilon(0.01));
344 1 : REQUIRE_EQ((c[2] - o[2]), Approx(0).epsilon(0.01));
345 1 : }
346 4 : }
347 :
348 20 : WHEN("testing geometry with 45/22.5 degree rotation and offset center of rotation")
349 20 : {
350 4 : real_t angle1 = pi_t / 4;
351 4 : real_t angle2 = pi_t / 2;
352 4 : RealVector_t offset(3);
353 4 : offset << 1, -2, -1;
354 :
355 4 : Geometry g(SourceToCenterOfRotation{s2c}, CenterOfRotationToDetector{c2d},
356 4 : VolumeData3D{Size3D{volCoeff}}, SinogramData3D{Size3D{detCoeff}},
357 4 : RotationAngles3D{Radian{angle1}, Radian{angle2}, Radian{0}},
358 4 : PrincipalPointOffset2D{0, 0}, RotationOffset3D{offset});
359 :
360 4 : THEN("a copy is the same")
361 4 : {
362 1 : Geometry gCopy(g);
363 1 : REQUIRE_EQ(gCopy, g);
364 1 : }
365 :
366 4 : THEN("then P and Pinv are inverse")
367 4 : {
368 1 : RealMatrix_t id(3, 3);
369 1 : id.setIdentity();
370 :
371 1 : auto P = g.getProjectionMatrix();
372 1 : auto Pinv = g.getInverseProjectionMatrix();
373 1 : RealMatrix_t result = (P * Pinv) - id;
374 1 : REQUIRE_EQ(result.sum(), Approx(0).epsilon(0.01));
375 1 : }
376 :
377 4 : THEN("then the rotation matrix is correct")
378 4 : {
379 1 : RealMatrix_t rot1(3, 3);
380 1 : rot1 << std::cos(angle1), 0, std::sin(angle1), 0, 1, 0, -std::sin(angle1), 0,
381 1 : std::cos(angle1);
382 1 : RealMatrix_t rot2(3, 3);
383 1 : rot2 << std::cos(angle2), -std::sin(angle2), 0, std::sin(angle2), std::cos(angle2),
384 1 : 0, 0, 0, 1;
385 1 : REQUIRE_EQ((g.getRotationMatrix() - rot1 * rot2).sum(), Approx(0));
386 1 : }
387 :
388 4 : THEN("the camera center is correct")
389 4 : {
390 1 : auto c = g.getCameraCenter();
391 1 : auto o = ddVol.getLocationOfOrigin();
392 :
393 1 : RealVector_t src(3);
394 1 : src[0] = 0;
395 1 : src[1] = 0;
396 1 : src[2] = -s2c;
397 1 : RealVector_t rotSrc = g.getRotationMatrix().transpose() * src + o + offset;
398 :
399 1 : REQUIRE_EQ((rotSrc - c).sum(), Approx(0).epsilon(0.01));
400 1 : }
401 4 : }
402 :
403 20 : WHEN("testing geometry with 45/22.5/12.25 degree rotation as a rotation matrix")
404 20 : {
405 4 : real_t angle1 = pi_t / 4;
406 4 : real_t angle2 = pi_t / 2;
407 4 : real_t angle3 = pi_t / 8;
408 4 : RealMatrix_t rot1(3, 3);
409 4 : rot1 << std::cos(angle1), 0, std::sin(angle1), 0, 1, 0, -std::sin(angle1), 0,
410 4 : std::cos(angle1);
411 4 : RealMatrix_t rot2(3, 3);
412 4 : rot2 << std::cos(angle2), -std::sin(angle2), 0, std::sin(angle2), std::cos(angle2), 0,
413 4 : 0, 0, 1;
414 4 : RealMatrix_t rot3(3, 3);
415 4 : rot3 << std::cos(angle3), 0, std::sin(angle3), 0, 1, 0, -std::sin(angle3), 0,
416 4 : std::cos(angle3);
417 4 : RealMatrix_t R = rot1 * rot2 * rot3;
418 :
419 4 : Geometry g(s2c, c2d, ddVol, ddDet, R);
420 :
421 4 : THEN("a copy is the same")
422 4 : {
423 1 : Geometry gCopy(g);
424 1 : REQUIRE_EQ(gCopy, g);
425 1 : }
426 :
427 4 : THEN("then P and Pinv are inverse")
428 4 : {
429 1 : RealMatrix_t id(3, 3);
430 1 : id.setIdentity();
431 :
432 1 : auto P = g.getProjectionMatrix();
433 1 : auto Pinv = g.getInverseProjectionMatrix();
434 1 : RealMatrix_t result = (P * Pinv) - id;
435 1 : REQUIRE_EQ(result.sum(), Approx(0).epsilon(0.01));
436 1 : }
437 :
438 4 : THEN("then the rotation matrix is correct")
439 4 : {
440 1 : REQUIRE_EQ((g.getRotationMatrix() - R).sum(), Approx(0));
441 1 : }
442 :
443 4 : THEN("the camera center is correct")
444 4 : {
445 1 : auto c = g.getCameraCenter();
446 1 : auto o = ddVol.getLocationOfOrigin();
447 :
448 1 : RealVector_t src(3);
449 1 : src[0] = 0;
450 1 : src[1] = 0;
451 1 : src[2] = -s2c;
452 1 : RealVector_t rotSrc = g.getRotationMatrix().transpose() * src + o;
453 :
454 1 : REQUIRE_EQ((rotSrc - c).sum(), Approx(0).epsilon(0.01));
455 1 : }
456 4 : }
457 20 : }
458 20 : }
459 :
460 : TEST_SUITE_END();
|