source: src/Fragmentation/Exporters/unittests/SphericalPointDistributionUnitTest.cpp@ 92a32a

Last change on this file since 92a32a was 92a32a, checked in by Frederik Heber <heber@…>, 9 years ago

Added "arbitrary" rotation for matchSphericalPointDistributionsTest_2().
  • Property mode set to 100644
File size: 16.3 KB
Line 
1/*
2 * Project: MoleCuilder
3 * Description: creates and alters molecular systems
4 * Copyright (C) 2014 Frederik Heber. All rights reserved.
5 *
6 *
7 * This file is part of MoleCuilder.
8 *
9 * MoleCuilder is free software: you can redistribute it and/or modify
10 * it under the terms of the GNU General Public License as published by
11 * the Free Software Foundation, either version 2 of the License, or
12 * (at your option) any later version.
13 *
14 * MoleCuilder 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. See the
17 * GNU General Public License for more details.
18 *
19 * You should have received a copy of the GNU General Public License
20 * along with MoleCuilder. If not, see <http://www.gnu.org/licenses/>.
21 */
22
23/*
24 * SphericalPointDistributionUnitTest.cpp
25 *
26 * Created on: May 29, 2014
27 * Author: heber
28 */
29
30// include config.h
31#ifdef HAVE_CONFIG_H
32#include <config.h>
33#endif
34
35using namespace std;
36
37#include <cppunit/CompilerOutputter.h>
38#include <cppunit/extensions/TestFactoryRegistry.h>
39#include <cppunit/ui/text/TestRunner.h>
40
41// include headers that implement a archive in simple text format
42#include <boost/archive/text_oarchive.hpp>
43#include <boost/archive/text_iarchive.hpp>
44
45#include "SphericalPointDistributionUnitTest.hpp"
46
47#include <boost/assign.hpp>
48#include <boost/math/quaternion.hpp>
49
50#include "CodePatterns/Assert.hpp"
51#include "CodePatterns/Log.hpp"
52
53#include "LinearAlgebra/Line.hpp"
54
55#include "Fragmentation/Exporters/SphericalPointDistribution.hpp"
56
57#include "LinearAlgebra/Line.hpp"
58
59#ifdef HAVE_TESTRUNNER
60#include "UnitTestMain.hpp"
61#endif /*HAVE_TESTRUNNER*/
62
63using namespace boost::assign;
64
65/********************************************** Test classes **************************************/
66
67// Registers the fixture into the 'registry'
68CPPUNIT_TEST_SUITE_REGISTRATION( SphericalPointDistributionTest );
69
70
71void SphericalPointDistributionTest::setUp()
72{
73 // failing asserts should be thrown
74 ASSERT_DO(Assert::Throw);
75
76 setVerbosity(5);
77}
78
79
80void SphericalPointDistributionTest::tearDown()
81{
82}
83
84void SphericalPointDistributionTest::QuaternionTest()
85{
86 Vector oldCenter(0.,1.,0.);
87 Vector newCenter(1.,0.,0.);
88
89 {
90 // setup quaternion
91 Vector RotationAxis = newCenter;
92 RotationAxis.VectorProduct(oldCenter);
93 RotationAxis.Normalize();
94 const double RotationAngle = oldCenter.Angle(newCenter)/(M_PI/2.);
95 // RotationAxis.Angle(oldCenter) - RotationAxis.Angle(newCenter);
96 boost::math::quaternion<double> q
97 (RotationAngle, RotationAxis[0], RotationAxis[1], RotationAxis[2]);
98 LOG(5, "DEBUG: RotationAxis is " << RotationAxis
99 << ", RotationAngle is " << RotationAngle);
100 LOG(5, "DEBUG: Quaternion describing rotation is " << q);
101 boost::math::quaternion<double> q_inverse =
102 boost::math::conj(q)/(boost::math::norm(q));
103 LOG(5, "DEBUG: Quaternion inverse is " << q_inverse);
104 boost::math::quaternion<double> identity(1,0,0,0);
105 const boost::math::quaternion<double> unity = q*q_inverse;
106 LOG(5, "DEBUG: q * q^-1 is " << unity);
107 CPPUNIT_ASSERT( boost::math::norm(unity - identity) < std::numeric_limits<double>::epsilon()*1e4);
108
109 // check that rotation works
110 boost::math::quaternion<double> p(0., newCenter[0], newCenter[1], newCenter[2]);
111 LOG(5, "DEBUG: Original newCenter is " << p);
112 p = p * q_inverse;
113 p = q * p;
114 LOG(5, "DEBUG: Rotated newCenter is " << p);
115 boost::math::quaternion<double> comparison(0., -oldCenter[0], oldCenter[1], oldCenter[2]);
116 LOG(5, "DEBUG: Difference norm is " << boost::math::norm(p - comparison));
117 CPPUNIT_ASSERT( boost::math::norm(p - comparison) < std::numeric_limits<double>::epsilon()*1e4);
118 }
119
120 // rotating with angle = 0 flips the vector unwantedly
121 {
122 // setup quaternion
123 Vector RotationAxis = newCenter;
124 RotationAxis.VectorProduct(oldCenter);
125 RotationAxis.Normalize();
126 const double RotationAngle = 0.;
127 // RotationAxis.Angle(oldCenter) - RotationAxis.Angle(newCenter);
128 boost::math::quaternion<double> q
129 (RotationAngle, RotationAxis[0], RotationAxis[1], RotationAxis[2]);
130 LOG(5, "DEBUG: RotationAxis is " << RotationAxis
131 << ", RotationAngle is " << RotationAngle);
132 LOG(5, "DEBUG: Quaternion describing rotation is " << q);
133 boost::math::quaternion<double> q_inverse =
134 boost::math::conj(q)/(boost::math::norm(q));
135 LOG(5, "DEBUG: Quaternion inverse is " << q_inverse);
136 boost::math::quaternion<double> identity(1,0,0,0);
137 const boost::math::quaternion<double> unity = q*q_inverse;
138 LOG(5, "DEBUG: q * q^-1 is " << unity);
139 CPPUNIT_ASSERT( boost::math::norm(unity - identity) < std::numeric_limits<double>::epsilon()*1e4);
140
141 // check that rotation works
142 boost::math::quaternion<double> p(0., newCenter[0], newCenter[1], newCenter[2]);
143 boost::math::quaternion<double> comparison(0., -newCenter[0], newCenter[1], newCenter[2]);
144 LOG(5, "DEBUG: Original newCenter is " << p);
145 p = p * q_inverse;
146 p = q * p;
147 LOG(5, "DEBUG: Rotated newCenter is " << p);
148 LOG(5, "DEBUG: Difference norm is " << boost::math::norm(p - comparison));
149 CPPUNIT_ASSERT( boost::math::norm(p - comparison) < std::numeric_limits<double>::epsilon()*1e4);
150 }
151}
152
153/** UnitTest for matchSphericalPointDistributions() with two points
154 */
155void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_2()
156{
157 SphericalPointDistribution SPD(1.);
158 // test with one point, matching trivially
159 {
160 SphericalPointDistribution::Polygon_t polygon;
161 polygon += Vector(1.,0.,0.);
162 SphericalPointDistribution::Polygon_t newpolygon =
163 SPD.get<2>();
164 SphericalPointDistribution::Polygon_t expected;
165 expected += Vector(-1.,0.,0.);
166 SphericalPointDistribution::Polygon_t remaining =
167 SphericalPointDistribution::matchSphericalPointDistributions(
168 polygon,
169 newpolygon);
170 CPPUNIT_ASSERT_EQUAL( expected, remaining );
171 }
172
173 // test with one point, just a flip of axis
174 {
175 SphericalPointDistribution::Polygon_t polygon;
176 polygon += Vector(0.,1.,0.);
177 SphericalPointDistribution::Polygon_t newpolygon =
178 SPD.get<2>();
179 SphericalPointDistribution::Polygon_t expected;
180 expected += Vector(0.,-1.,0.);
181 SphericalPointDistribution::Polygon_t remaining =
182 SphericalPointDistribution::matchSphericalPointDistributions(
183 polygon,
184 newpolygon);
185 CPPUNIT_ASSERT_EQUAL( expected, remaining );
186 }
187
188 // test with one point, just a flip to another axis
189 {
190 SphericalPointDistribution::Polygon_t polygon;
191 polygon += Vector(0.,0.,-1.);
192 SphericalPointDistribution::Polygon_t newpolygon =
193 SPD.get<2>();
194 SphericalPointDistribution::Polygon_t expected;
195 expected += Vector(0.,0.,1.);
196 SphericalPointDistribution::Polygon_t remaining =
197 SphericalPointDistribution::matchSphericalPointDistributions(
198 polygon,
199 newpolygon);
200 CPPUNIT_ASSERT_EQUAL( expected, remaining );
201 }
202
203 // test with one point, full rotation
204 {
205 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
206 SphericalPointDistribution::Polygon_t polygon;
207 polygon += RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI);
208 SphericalPointDistribution::Polygon_t newpolygon =
209 SPD.get<2>();
210 SphericalPointDistribution::Polygon_t expected;
211 expected += RotationAxis.rotateVector(Vector(-1.,0.,0.), 47.6/180*M_PI);
212 SphericalPointDistribution::Polygon_t remaining =
213 SphericalPointDistribution::matchSphericalPointDistributions(
214 polygon,
215 newpolygon);
216 CPPUNIT_ASSERT_EQUAL( expected, remaining );
217 }
218}
219
220/** UnitTest for matchSphericalPointDistributions() with three points
221 */
222void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_3()
223{
224 SphericalPointDistribution SPD(1.);
225
226 // test with one point, matching trivially
227 {
228 SphericalPointDistribution::Polygon_t polygon;
229 polygon += Vector(1.,0.,0.);
230 SphericalPointDistribution::Polygon_t newpolygon =
231 SPD.get<3>();
232 SphericalPointDistribution::Polygon_t expected = newpolygon;
233 expected.pop_front(); // remove first point
234 SphericalPointDistribution::Polygon_t remaining =
235 SphericalPointDistribution::matchSphericalPointDistributions(
236 polygon,
237 newpolygon);
238 CPPUNIT_ASSERT_EQUAL( expected, remaining );
239 }
240
241 // test with one point, just a flip of x and y axis
242 {
243 SphericalPointDistribution::Polygon_t polygon;
244 polygon += Vector(0.,1.,0.);
245 SphericalPointDistribution::Polygon_t newpolygon =
246 SPD.get<3>();
247 SphericalPointDistribution::Polygon_t expected = newpolygon;
248 expected.pop_front(); // remove first point
249 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
250 iter != expected.end(); ++iter) {
251 std::swap((*iter)[0], (*iter)[1]);
252 (*iter)[0] *= -1.;
253 }
254 SphericalPointDistribution::Polygon_t remaining =
255 SphericalPointDistribution::matchSphericalPointDistributions(
256 polygon,
257 newpolygon);
258 CPPUNIT_ASSERT_EQUAL( expected, remaining );
259 }
260}
261
262/** UnitTest for matchSphericalPointDistributions() with four points
263 */
264void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_4()
265{
266 SphericalPointDistribution SPD(1.);
267
268 // test with one point, matching trivially
269 {
270 SphericalPointDistribution::Polygon_t polygon;
271 polygon += Vector(1.,0.,0.);
272 SphericalPointDistribution::Polygon_t newpolygon =
273 SPD.get<4>();
274 SphericalPointDistribution::Polygon_t expected = newpolygon;
275 expected.pop_front(); // remove first point
276 SphericalPointDistribution::Polygon_t remaining =
277 SphericalPointDistribution::matchSphericalPointDistributions(
278 polygon,
279 newpolygon);
280 CPPUNIT_ASSERT_EQUAL( expected, remaining );
281 }
282
283 // test with one point, just a flip of axis
284 {
285 SphericalPointDistribution::Polygon_t polygon;
286 polygon += Vector(0.,1.,0.);
287 SphericalPointDistribution::Polygon_t newpolygon =
288 SPD.get<4>();
289 SphericalPointDistribution::Polygon_t expected = newpolygon;
290 expected.pop_front(); // remove first point
291 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
292 iter != expected.end(); ++iter) {
293 std::swap((*iter)[0], (*iter)[1]);
294 (*iter)[0] *= -1.;
295 }
296 SphericalPointDistribution::Polygon_t remaining =
297 SphericalPointDistribution::matchSphericalPointDistributions(
298 polygon,
299 newpolygon);
300 CPPUNIT_ASSERT_EQUAL( expected, remaining );
301 }
302}
303
304/** UnitTest for matchSphericalPointDistributions() with five points
305 */
306void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_5()
307{
308 SphericalPointDistribution SPD(1.);
309
310 // test with one point, matching trivially
311 {
312 SphericalPointDistribution::Polygon_t polygon;
313 polygon += Vector(1.,0.,0.);
314 SphericalPointDistribution::Polygon_t newpolygon =
315 SPD.get<5>();
316 SphericalPointDistribution::Polygon_t expected = newpolygon;
317 expected.pop_front(); // remove first point
318 SphericalPointDistribution::Polygon_t remaining =
319 SphericalPointDistribution::matchSphericalPointDistributions(
320 polygon,
321 newpolygon);
322 CPPUNIT_ASSERT_EQUAL( expected, remaining );
323 }
324
325 // test with one point, just a flip of axis
326 {
327 SphericalPointDistribution::Polygon_t polygon;
328 polygon += Vector(0.,1.,0.);
329 SphericalPointDistribution::Polygon_t newpolygon =
330 SPD.get<5>();
331 SphericalPointDistribution::Polygon_t expected = newpolygon;
332 expected.pop_front(); // remove first point
333 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
334 iter != expected.end(); ++iter) {
335 std::swap((*iter)[0], (*iter)[1]);
336 (*iter)[0] *= -1.;
337 }
338 SphericalPointDistribution::Polygon_t remaining =
339 SphericalPointDistribution::matchSphericalPointDistributions(
340 polygon,
341 newpolygon);
342 CPPUNIT_ASSERT_EQUAL( expected, remaining );
343 }
344}
345
346/** UnitTest for matchSphericalPointDistributions() with six points
347 */
348void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_6()
349{
350 SphericalPointDistribution SPD(1.);
351
352 // test with one point, matching trivially
353 {
354 SphericalPointDistribution::Polygon_t polygon;
355 polygon += Vector(1.,0.,0.);
356 SphericalPointDistribution::Polygon_t newpolygon =
357 SPD.get<6>();
358 SphericalPointDistribution::Polygon_t expected = newpolygon;
359 expected.pop_front(); // remove first point
360 SphericalPointDistribution::Polygon_t remaining =
361 SphericalPointDistribution::matchSphericalPointDistributions(
362 polygon,
363 newpolygon);
364 CPPUNIT_ASSERT_EQUAL( expected, remaining );
365 }
366
367 // test with one point, just a flip of axis
368 {
369 SphericalPointDistribution::Polygon_t polygon;
370 polygon += Vector(0.,1.,0.);
371 SphericalPointDistribution::Polygon_t newpolygon =
372 SPD.get<6>();
373 SphericalPointDistribution::Polygon_t expected = newpolygon;
374 expected.pop_front(); // remove first point
375 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
376 iter != expected.end(); ++iter) {
377 std::swap((*iter)[0], (*iter)[1]);
378 (*iter)[0] *= -1.;
379 }
380 SphericalPointDistribution::Polygon_t remaining =
381 SphericalPointDistribution::matchSphericalPointDistributions(
382 polygon,
383 newpolygon);
384 CPPUNIT_ASSERT_EQUAL( expected, remaining );
385 }
386}
387
388/** UnitTest for matchSphericalPointDistributions() with seven points
389 */
390void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_7()
391{
392 SphericalPointDistribution SPD(1.);
393
394 // test with one point, matching trivially
395 {
396 SphericalPointDistribution::Polygon_t polygon;
397 polygon += Vector(1.,0.,0.);
398 SphericalPointDistribution::Polygon_t newpolygon =
399 SPD.get<7>();
400 SphericalPointDistribution::Polygon_t expected = newpolygon;
401 expected.pop_front(); // remove first point
402 SphericalPointDistribution::Polygon_t remaining =
403 SphericalPointDistribution::matchSphericalPointDistributions(
404 polygon,
405 newpolygon);
406 CPPUNIT_ASSERT_EQUAL( expected, remaining );
407 }
408
409 // test with one point, just a flip of axis
410 {
411 SphericalPointDistribution::Polygon_t polygon;
412 polygon += Vector(0.,1.,0.);
413 SphericalPointDistribution::Polygon_t newpolygon =
414 SPD.get<7>();
415 SphericalPointDistribution::Polygon_t expected = newpolygon;
416 expected.pop_front(); // remove first point
417 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
418 iter != expected.end(); ++iter) {
419 std::swap((*iter)[0], (*iter)[1]);
420 (*iter)[0] *= -1.;
421 }
422 SphericalPointDistribution::Polygon_t remaining =
423 SphericalPointDistribution::matchSphericalPointDistributions(
424 polygon,
425 newpolygon);
426 CPPUNIT_ASSERT_EQUAL( expected, remaining );
427 }
428}
429
430/** UnitTest for matchSphericalPointDistributions() with eight points
431 */
432void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_8()
433{
434 SphericalPointDistribution SPD(1.);
435
436 // test with one point, matching trivially
437 {
438 SphericalPointDistribution::Polygon_t polygon;
439 polygon += Vector(1.,0.,0.);
440 SphericalPointDistribution::Polygon_t newpolygon =
441 SPD.get<8>();
442 SphericalPointDistribution::Polygon_t expected = newpolygon;
443 expected.pop_front(); // remove first point
444 SphericalPointDistribution::Polygon_t remaining =
445 SphericalPointDistribution::matchSphericalPointDistributions(
446 polygon,
447 newpolygon);
448 CPPUNIT_ASSERT_EQUAL( expected, remaining );
449 }
450
451 // test with one point, just a flip of axis
452 {
453 SphericalPointDistribution::Polygon_t polygon;
454 polygon += Vector(0.,1.,0.);
455 SphericalPointDistribution::Polygon_t newpolygon =
456 SPD.get<8>();
457 SphericalPointDistribution::Polygon_t expected = newpolygon;
458 expected.pop_front(); // remove first point
459 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
460 iter != expected.end(); ++iter) {
461 std::swap((*iter)[0], (*iter)[1]);
462 (*iter)[0] *= -1.;
463 }
464 SphericalPointDistribution::Polygon_t remaining =
465 SphericalPointDistribution::matchSphericalPointDistributions(
466 polygon,
467 newpolygon);
468 CPPUNIT_ASSERT_EQUAL( expected, remaining );
469 }
470}
Note: See TracBrowser for help on using the repository browser.