source: src/Fragmentation/Exporters/unittests/SphericalPointDistributionUnitTest.cpp@ 0983e6

Action_Thermostats Add_AtomRandomPerturbation Add_FitFragmentPartialChargesAction Add_RotateAroundBondAction Add_SelectAtomByNameAction Adding_Graph_to_ChangeBondActions Adding_MD_integration_tests Adding_StructOpt_integration_tests Automaking_mpqc_open AutomationFragmentation_failures Candidate_v1.5.4 Candidate_v1.6.0 Candidate_v1.6.1 ChangeBugEmailaddress ChangingTestPorts ChemicalSpaceEvaluator Combining_Subpackages Debian_Package_split Debian_package_split_molecuildergui_only Disabling_MemDebug Docu_Python_wait EmpiricalPotential_contain_HomologyGraph EmpiricalPotential_contain_HomologyGraph_documentation Enable_parallel_make_install Enhance_userguide Enhanced_StructuralOptimization Enhanced_StructuralOptimization_continued Example_ManyWaysToTranslateAtom Exclude_Hydrogens_annealWithBondGraph FitPartialCharges_GlobalError Fix_ChargeSampling_PBC Fix_ChronosMutex Fix_FitPartialCharges Fix_FitPotential_needs_atomicnumbers Fix_ForceAnnealing Fix_IndependentFragmentGrids Fix_ParseParticles Fix_ParseParticles_split_forward_backward_Actions Fix_StatusMsg Fix_StepWorldTime_single_argument Fix_Verbose_Codepatterns ForceAnnealing_goodresults ForceAnnealing_oldresults ForceAnnealing_tocheck ForceAnnealing_with_BondGraph ForceAnnealing_with_BondGraph_continued ForceAnnealing_with_BondGraph_continued_betteresults ForceAnnealing_with_BondGraph_contraction-expansion GeometryObjects Gui_displays_atomic_force_velocity IndependentFragmentGrids IndependentFragmentGrids_IndividualZeroInstances IndependentFragmentGrids_IntegrationTest IndependentFragmentGrids_Sole_NN_Calculation JobMarket_RobustOnKillsSegFaults JobMarket_StableWorkerPool JobMarket_unresolvable_hostname_fix ODR_violation_mpqc_open PartialCharges_OrthogonalSummation PythonUI_with_named_parameters QtGui_reactivate_TimeChanged_changes Recreated_GuiChecks RotateToPrincipalAxisSystem_UndoRedo SaturateAtoms_findBestMatching StoppableMakroAction Subpackage_CodePatterns Subpackage_JobMarket Subpackage_LinearAlgebra Subpackage_levmar Subpackage_mpqc_open Subpackage_vmg ThirdParty_MPQC_rebuilt_buildsystem TrajectoryDependenant_MaxOrder TremoloParser_IncreasedPrecision TremoloParser_MultipleTimesteps Ubuntu_1604_changes stable
Last change on this file since 0983e6 was 653cea, checked in by Frederik Heber <heber@…>, 8 years ago

SphericalPointDistribution is now working with bond degree weights.

  • recurseMatching() now works on IndexTupleList_t.
  • also rewrote calculatePairwiseDistances() and calculateErrorOfMatching().
  • L1THRESHOLD in recurseMatching() moved over to class body.
  • increased verbosity level of ...Matching() functions by one, added note on eventually chosen matching and why.
  • we assert that bestL2 is not too large.
  • FIX: calculateErrorOfMatching() did not use absolute value of gap for L1 error.
  • TESTFIX: Using limited accuracy on point coordinates.
  • TESTS: Regresssion test FragmentMolecule-cycles working again.
  • Property mode set to 100644
File size: 53.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/** due to root-taking in function we only have limited numerical precision,
71 * basically half of the double range.
72 */
73const double CenterAccuracy = sqrt(std::numeric_limits<double>::epsilon()*1e2);
74
75void SphericalPointDistributionTest::setUp()
76{
77 // failing asserts should be thrown
78 ASSERT_DO(Assert::Throw);
79
80 setVerbosity(6);
81}
82
83
84void SphericalPointDistributionTest::tearDown()
85{
86}
87
88
89/** UnitTest for matchSphericalPointDistributions() with two points
90 */
91void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_2()
92{
93 SphericalPointDistribution SPD(1.);
94 // test with one point, matching trivially
95 {
96 SphericalPointDistribution::WeightedPolygon_t polygon;
97 polygon += std::make_pair(Vector(1.,0.,0.), 1);
98 SphericalPointDistribution::Polygon_t newpolygon =
99 SPD.get<2>();
100 SphericalPointDistribution::Polygon_t expected;
101 expected += Vector(-1.,0.,0.);
102 SphericalPointDistribution::Polygon_t remaining =
103 SphericalPointDistribution::matchSphericalPointDistributions(
104 polygon,
105 newpolygon);
106 CPPUNIT_ASSERT_EQUAL( expected, remaining );
107 }
108
109 // test with one point, just a flip of axis
110 {
111 SphericalPointDistribution::WeightedPolygon_t polygon;
112 polygon += std::make_pair( Vector(0.,1.,0.), 1);
113 SphericalPointDistribution::Polygon_t newpolygon =
114 SPD.get<2>();
115 SphericalPointDistribution::Polygon_t expected;
116 expected += Vector(0.,-1.,0.);
117 SphericalPointDistribution::Polygon_t remaining =
118 SphericalPointDistribution::matchSphericalPointDistributions(
119 polygon,
120 newpolygon);
121 CPPUNIT_ASSERT_EQUAL( expected, remaining );
122 }
123
124 // test with one point, just a flip to another axis
125 {
126 SphericalPointDistribution::WeightedPolygon_t polygon;
127 polygon += std::make_pair( Vector(0.,0.,-1.), 1);
128 SphericalPointDistribution::Polygon_t newpolygon =
129 SPD.get<2>();
130 SphericalPointDistribution::Polygon_t expected;
131 expected += Vector(0.,0.,1.);
132 SphericalPointDistribution::Polygon_t remaining =
133 SphericalPointDistribution::matchSphericalPointDistributions(
134 polygon,
135 newpolygon);
136 CPPUNIT_ASSERT_EQUAL( expected, remaining );
137 }
138
139 // test with one point, full rotation
140 {
141 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
142 SphericalPointDistribution::WeightedPolygon_t polygon;
143 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
144 SphericalPointDistribution::Polygon_t newpolygon =
145 SPD.get<2>();
146 SphericalPointDistribution::Polygon_t expected;
147 expected += RotationAxis.rotateVector(Vector(-1.,0.,0.), 47.6/180*M_PI);
148 SphericalPointDistribution::Polygon_t remaining =
149 SphericalPointDistribution::matchSphericalPointDistributions(
150 polygon,
151 newpolygon);
152 CPPUNIT_ASSERT_EQUAL( expected, remaining );
153 }
154}
155
156void perturbPolygon(
157 SphericalPointDistribution::WeightedPolygon_t &_polygon,
158 double _amplitude
159 )
160{
161 for (SphericalPointDistribution::WeightedPolygon_t::iterator iter = _polygon.begin();
162 iter != _polygon.end(); ++iter) {
163 Vector perturber;
164 perturber.GetOneNormalVector(iter->first);
165 perturber.Scale(_amplitude);
166 iter->first = iter->first + perturber;
167 (iter->first).Normalize();
168 }
169}
170
171static
172bool areEqualToWithinBounds(
173 const SphericalPointDistribution::Polygon_t &_polygon,
174 const SphericalPointDistribution::Polygon_t &_otherpolygon,
175 double _amplitude
176 )
177{
178 // same size?
179 if (_polygon.size() != _otherpolygon.size())
180 return false;
181 // same points ? We just check witrh trivial mapping, nothing fancy ...
182 bool status = true;
183 SphericalPointDistribution::Polygon_t::const_iterator iter = _polygon.begin();
184 SphericalPointDistribution::Polygon_t::const_iterator otheriter = _otherpolygon.begin();
185 for (; iter != _polygon.end(); ++iter, ++otheriter) {
186 status &= (*iter - *otheriter).Norm() < _amplitude;
187 }
188 return status;
189}
190
191/** UnitTest for areEqualToWithinBounds()
192 */
193void SphericalPointDistributionTest::areEqualToWithinBoundsTest()
194{
195 // test with no points
196 {
197 SphericalPointDistribution::Polygon_t polygon;
198 SphericalPointDistribution::Polygon_t expected = polygon;
199 CPPUNIT_ASSERT( areEqualToWithinBounds(polygon, expected, std::numeric_limits<double>::epsilon()*1e2) );
200 }
201 // test with one point
202 {
203 SphericalPointDistribution::Polygon_t polygon;
204 polygon += Vector(1.,0.,0.);
205 SphericalPointDistribution::Polygon_t expected = polygon;
206 CPPUNIT_ASSERT( areEqualToWithinBounds(polygon, expected, std::numeric_limits<double>::epsilon()*1e2) );
207 }
208 // test with two points
209 {
210 SphericalPointDistribution::Polygon_t polygon;
211 polygon += Vector(1.,0.,0.);
212 polygon += Vector(0.,1.,0.);
213 SphericalPointDistribution::Polygon_t expected = polygon;
214 CPPUNIT_ASSERT( areEqualToWithinBounds(polygon, expected, std::numeric_limits<double>::epsilon()*1e2) );
215 }
216
217 // test with two points in different order: THIS GOES WRONG: We only check trivially
218 {
219 SphericalPointDistribution::Polygon_t polygon;
220 polygon += Vector(1.,0.,0.);
221 polygon += Vector(0.,1.,0.);
222 SphericalPointDistribution::Polygon_t expected;
223 expected += Vector(0.,1.,0.);
224 expected += Vector(1.,0.,0.);
225 CPPUNIT_ASSERT( !areEqualToWithinBounds(polygon, expected, std::numeric_limits<double>::epsilon()*1e2) );
226 }
227
228 // test with two different points
229 {
230 SphericalPointDistribution::Polygon_t polygon;
231 polygon += Vector(1.,0.,0.);
232 polygon += Vector(0.,1.,0.);
233 SphericalPointDistribution::Polygon_t expected;
234 expected += Vector(1.01,0.,0.);
235 expected += Vector(0.,1.,0.);
236 CPPUNIT_ASSERT( areEqualToWithinBounds(polygon, expected, 0.05) );
237 CPPUNIT_ASSERT( !areEqualToWithinBounds(polygon, expected, 0.005) );
238 }
239
240 // test with different number of points
241 {
242 SphericalPointDistribution::Polygon_t polygon;
243 polygon += Vector(1.,0.,0.);
244 polygon += Vector(0.,1.,0.);
245 SphericalPointDistribution::Polygon_t expected;
246 expected += Vector(0.,1.,0.);
247 CPPUNIT_ASSERT( !areEqualToWithinBounds(polygon, expected, 0.05) );
248 }
249}
250
251/** UnitTest for joinPoints()
252 */
253void SphericalPointDistributionTest::joinPointsTest()
254{
255 // test with simple configuration of three points
256 {
257 SphericalPointDistribution::Polygon_t newpolygon;
258 newpolygon += Vector(1.,0.,0.);
259 newpolygon += Vector(0.,1.,0.);
260 newpolygon += Vector(0.,0.,1.);
261 SphericalPointDistribution::Polygon_t expectedpolygon = newpolygon;
262 SphericalPointDistribution::IndexTupleList_t matching;
263 matching += SphericalPointDistribution::IndexList_t(1,0);
264 matching += SphericalPointDistribution::IndexList_t(1,1);
265 matching += SphericalPointDistribution::IndexList_t(1,2);
266 SphericalPointDistribution::IndexList_t IndexList =
267 SphericalPointDistribution::joinPoints(
268 newpolygon,
269 SphericalPointDistribution::VectorArray_t(newpolygon.begin(), newpolygon.end()),
270 matching);
271 SphericalPointDistribution::IndexList_t expected;
272 expected += 0,1,2;
273 CPPUNIT_ASSERT_EQUAL( expected, IndexList );
274 CPPUNIT_ASSERT_EQUAL( expectedpolygon, newpolygon );
275 }
276
277 // test with simple configuration of three points, only two are picked
278 {
279 SphericalPointDistribution::Polygon_t newpolygon;
280 newpolygon += Vector(1.,0.,0.);
281 newpolygon += Vector(0.,1.,0.);
282 newpolygon += Vector(0.,0.,1.);
283 SphericalPointDistribution::Polygon_t expectedpolygon = newpolygon;
284 SphericalPointDistribution::IndexTupleList_t matching;
285 matching += SphericalPointDistribution::IndexList_t(1,1);
286 matching += SphericalPointDistribution::IndexList_t(1,2);
287 SphericalPointDistribution::IndexList_t IndexList =
288 SphericalPointDistribution::joinPoints(
289 newpolygon,
290 SphericalPointDistribution::VectorArray_t(newpolygon.begin(), newpolygon.end()),
291 matching);
292 SphericalPointDistribution::IndexList_t expected;
293 expected += 1,2;
294 CPPUNIT_ASSERT_EQUAL( expected, IndexList );
295 CPPUNIT_ASSERT_EQUAL( expectedpolygon, newpolygon );
296 }
297
298 // test with simple configuration of three points, two are joined
299 {
300 SphericalPointDistribution::Polygon_t newpolygon;
301 newpolygon += Vector(1.,0.,0.);
302 newpolygon += Vector(0.,1.,0.);
303 newpolygon += Vector(0.,0.,1.);
304 SphericalPointDistribution::Polygon_t expectedpolygon;
305 expectedpolygon += Vector(1.,0.,0.);
306 expectedpolygon += Vector(0.,M_SQRT1_2,M_SQRT1_2);
307 SphericalPointDistribution::IndexTupleList_t matching;
308 SphericalPointDistribution::IndexList_t joined;
309 joined += 1,2;
310 matching += SphericalPointDistribution::IndexList_t(1,0);
311 matching += joined;
312 SphericalPointDistribution::IndexList_t IndexList =
313 SphericalPointDistribution::joinPoints(
314 newpolygon,
315 SphericalPointDistribution::VectorArray_t(newpolygon.begin(), newpolygon.end()),
316 matching);
317 SphericalPointDistribution::IndexList_t expected;
318 expected += 0,1;
319 CPPUNIT_ASSERT_EQUAL( expected, IndexList );
320 CPPUNIT_ASSERT_EQUAL( expectedpolygon, newpolygon );
321 }
322
323 // test with simple configuration of six points, two are joined, jumbled indices
324 {
325 SphericalPointDistribution::Polygon_t newpolygon;
326 newpolygon += Vector(1.,0.,1.);
327 newpolygon += Vector(1.,0.,0.);
328 newpolygon += Vector(1.,1.,0.);
329 newpolygon += Vector(0.,1.,0.);
330 newpolygon += Vector(0.,0.,1.);
331 newpolygon += Vector(1.,0.,1.);
332 SphericalPointDistribution::Polygon_t expectedpolygon;
333 expectedpolygon += Vector(1.,0.,1.);
334 expectedpolygon += Vector(1.,0.,0.);
335 expectedpolygon += Vector(1.,1.,0.);
336 expectedpolygon += Vector(1.,0.,1.);
337 expectedpolygon += Vector(0.,M_SQRT1_2,M_SQRT1_2); // new centers go last
338 SphericalPointDistribution::IndexTupleList_t matching;
339 SphericalPointDistribution::IndexList_t joined;
340 joined += 3,4;
341 matching += SphericalPointDistribution::IndexList_t(1,1);
342 matching += joined;
343 SphericalPointDistribution::IndexList_t IndexList =
344 SphericalPointDistribution::joinPoints(
345 newpolygon,
346 SphericalPointDistribution::VectorArray_t(newpolygon.begin(), newpolygon.end()),
347 matching);
348 SphericalPointDistribution::IndexList_t expected;
349 expected += 1,4;
350 CPPUNIT_ASSERT_EQUAL( expected, IndexList );
351 CPPUNIT_ASSERT_EQUAL( expectedpolygon, newpolygon );
352 }
353}
354
355/** UnitTest for matchSphericalPointDistributions() with three points
356 */
357void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_3()
358{
359 SphericalPointDistribution SPD(1.);
360
361 // test with one point, matching trivially
362 {
363 SphericalPointDistribution::WeightedPolygon_t polygon;
364 polygon += std::make_pair( Vector(1.,0.,0.), 1);
365 SphericalPointDistribution::Polygon_t newpolygon =
366 SPD.get<3>();
367 SphericalPointDistribution::Polygon_t expected = newpolygon;
368 expected.pop_front(); // remove first point
369 SphericalPointDistribution::Polygon_t remaining =
370 SphericalPointDistribution::matchSphericalPointDistributions(
371 polygon,
372 newpolygon);
373 CPPUNIT_ASSERT_EQUAL( expected, remaining );
374 }
375
376 // test with one point, just a flip of x and y axis
377 {
378 SphericalPointDistribution::WeightedPolygon_t polygon;
379 polygon += std::make_pair( Vector(0.,1.,0.), 1);
380 SphericalPointDistribution::Polygon_t newpolygon =
381 SPD.get<3>();
382 SphericalPointDistribution::Polygon_t expected = newpolygon;
383 expected.pop_front(); // remove first point
384 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
385 iter != expected.end(); ++iter) {
386 std::swap((*iter)[0], (*iter)[1]);
387 (*iter)[0] *= -1.;
388 }
389 SphericalPointDistribution::Polygon_t remaining =
390 SphericalPointDistribution::matchSphericalPointDistributions(
391 polygon,
392 newpolygon);
393 CPPUNIT_ASSERT_EQUAL( expected, remaining );
394 }
395
396 // test with two points, matching trivially
397 {
398 SphericalPointDistribution::WeightedPolygon_t polygon;
399 polygon += std::make_pair( Vector(1.,0.,0.), 1);
400 polygon += std::make_pair( Vector(-0.5, sqrt(3)*0.5,0.), 1);
401 SphericalPointDistribution::Polygon_t newpolygon =
402 SPD.get<3>();
403 SphericalPointDistribution::Polygon_t expected = newpolygon;
404 expected.pop_front(); // remove first point
405 expected.pop_front(); // remove second point
406 SphericalPointDistribution::Polygon_t remaining =
407 SphericalPointDistribution::matchSphericalPointDistributions(
408 polygon,
409 newpolygon);
410 CPPUNIT_ASSERT_EQUAL( expected, remaining );
411 // also slightly perturbed
412 const double amplitude = 0.05;
413 perturbPolygon(polygon, amplitude);
414 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
415 }
416
417 // test with two points, full rotation
418 {
419 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
420 SphericalPointDistribution::WeightedPolygon_t polygon;
421 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
422 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-0.5, sqrt(3)*0.5,0.), 47.6/180*M_PI), 1);
423 SphericalPointDistribution::Polygon_t newpolygon =
424 SPD.get<3>();
425 SphericalPointDistribution::Polygon_t expected = newpolygon;
426 expected.pop_front(); // remove first point
427 expected.pop_front(); // remove second point
428 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
429 iter != expected.end(); ++iter)
430 *iter = RotationAxis.rotateVector(*iter, 47.6/180*M_PI);
431 SphericalPointDistribution::Polygon_t remaining =
432 SphericalPointDistribution::matchSphericalPointDistributions(
433 polygon,
434 newpolygon);
435 CPPUNIT_ASSERT_EQUAL( expected, remaining );
436 // also slightly perturbed
437 const double amplitude = 0.05;
438 perturbPolygon(polygon, amplitude);
439 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
440 }
441
442 // test with three points, matching trivially
443 {
444 SphericalPointDistribution::WeightedPolygon_t polygon;
445 polygon += std::make_pair( Vector(1.,0.,0.), 1);
446 polygon += std::make_pair( Vector(-0.5, sqrt(3)*0.5,0.), 1);
447 polygon += std::make_pair( Vector(-0.5, -sqrt(3)*0.5,0.), 1);
448 SphericalPointDistribution::Polygon_t newpolygon =
449 SPD.get<3>();
450 SphericalPointDistribution::Polygon_t expected; // empty cause none are vacant
451 SphericalPointDistribution::Polygon_t remaining =
452 SphericalPointDistribution::matchSphericalPointDistributions(
453 polygon,
454 newpolygon);
455 CPPUNIT_ASSERT_EQUAL( expected, remaining );
456 // also slightly perturbed
457 const double amplitude = 0.05;
458 perturbPolygon(polygon, amplitude);
459 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
460 }
461
462
463 // test with three points, full rotation
464 {
465 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
466 SphericalPointDistribution::WeightedPolygon_t polygon;
467 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
468 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-0.5, sqrt(3)*0.5,0.), 47.6/180*M_PI), 1);
469 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-0.5, -sqrt(3)*0.5,0.), 47.6/180*M_PI), 1);
470 SphericalPointDistribution::Polygon_t newpolygon =
471 SPD.get<3>();
472 SphericalPointDistribution::Polygon_t expected; // empty cause none are vacant
473 SphericalPointDistribution::Polygon_t remaining =
474 SphericalPointDistribution::matchSphericalPointDistributions(
475 polygon,
476 newpolygon);
477 CPPUNIT_ASSERT_EQUAL( expected, remaining );
478 // also slightly perturbed
479 const double amplitude = 0.05;
480 perturbPolygon(polygon, amplitude);
481 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
482 }
483}
484
485/** UnitTest for matchSphericalPointDistributions() with four points
486 */
487void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_4()
488{
489 SphericalPointDistribution SPD(1.);
490
491 // test with one point, matching trivially
492 {
493 SphericalPointDistribution::WeightedPolygon_t polygon;
494 polygon += std::make_pair( Vector(1.,0.,0.), 1);
495 SphericalPointDistribution::Polygon_t newpolygon =
496 SPD.get<4>();
497 SphericalPointDistribution::Polygon_t expected = newpolygon;
498 expected.pop_front(); // remove first point
499 SphericalPointDistribution::Polygon_t remaining =
500 SphericalPointDistribution::matchSphericalPointDistributions(
501 polygon,
502 newpolygon);
503 CPPUNIT_ASSERT_EQUAL( expected, remaining );
504 }
505
506 // test with one point, just a flip of axis
507 {
508 SphericalPointDistribution::WeightedPolygon_t polygon;
509 polygon += std::make_pair( Vector(0.,1.,0.), 1);
510 SphericalPointDistribution::Polygon_t newpolygon =
511 SPD.get<4>();
512 SphericalPointDistribution::Polygon_t expected = newpolygon;
513 expected.pop_front(); // remove first point
514 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
515 iter != expected.end(); ++iter) {
516 std::swap((*iter)[0], (*iter)[1]);
517 (*iter)[0] *= -1.;
518 }
519 SphericalPointDistribution::Polygon_t remaining =
520 SphericalPointDistribution::matchSphericalPointDistributions(
521 polygon,
522 newpolygon);
523 CPPUNIT_ASSERT_EQUAL( expected, remaining );
524 }
525
526 // test with two points, matching trivially
527 {
528 SphericalPointDistribution::WeightedPolygon_t polygon;
529 polygon += std::make_pair( Vector(1.,0.,0.), 1);
530 polygon += std::make_pair( Vector(-1./3.0, 2.0*M_SQRT2/3.0,0.), 1);
531 SphericalPointDistribution::Polygon_t newpolygon =
532 SPD.get<4>();
533 SphericalPointDistribution::Polygon_t expected = newpolygon;
534 expected.pop_front(); // remove first point
535 expected.pop_front(); // remove second point
536 SphericalPointDistribution::Polygon_t remaining =
537 SphericalPointDistribution::matchSphericalPointDistributions(
538 polygon,
539 newpolygon);
540 CPPUNIT_ASSERT_EQUAL( expected, remaining );
541 // also slightly perturbed
542 const double amplitude = 0.05;
543 perturbPolygon(polygon, amplitude);
544 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
545 }
546
547 // test with two points, matching trivially, also with slightly perturbed
548 {
549 SphericalPointDistribution::WeightedPolygon_t polygon;
550 polygon += std::make_pair( Vector(1.,0.,0.), 1);
551 polygon += std::make_pair( Vector(-1./3.0, 2.0*M_SQRT2/3.0,0.), 1);
552 SphericalPointDistribution::Polygon_t newpolygon =
553 SPD.get<4>();
554 SphericalPointDistribution::Polygon_t expected = newpolygon;
555 expected.pop_front(); // remove first point
556 expected.pop_front(); // remove second point
557 SphericalPointDistribution::Polygon_t remaining =
558 SphericalPointDistribution::matchSphericalPointDistributions(
559 polygon,
560 newpolygon);
561 CPPUNIT_ASSERT_EQUAL( expected, remaining );
562 // also slightly perturbed
563 const double amplitude = 0.05;
564 perturbPolygon(polygon, amplitude);
565 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
566 }
567
568 // test with two points, full rotation
569 {
570 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
571 SphericalPointDistribution::WeightedPolygon_t polygon;
572 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
573 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1./3.0, 2.0*M_SQRT2/3.0,0.), 47.6/180*M_PI), 1);
574 SphericalPointDistribution::Polygon_t newpolygon =
575 SPD.get<4>();
576 SphericalPointDistribution::Polygon_t expected = newpolygon;
577 expected.pop_front(); // remove first point
578 expected.pop_front(); // remove second point
579 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
580 iter != expected.end(); ++iter)
581 *iter = RotationAxis.rotateVector(*iter, 47.6/180*M_PI);
582 SphericalPointDistribution::Polygon_t remaining =
583 SphericalPointDistribution::matchSphericalPointDistributions(
584 polygon,
585 newpolygon);
586 CPPUNIT_ASSERT_EQUAL( expected, remaining );
587 // also slightly perturbed
588 const double amplitude = 0.05;
589 perturbPolygon(polygon, amplitude);
590 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
591 }
592
593 // test with three points, matching trivially
594 {
595 SphericalPointDistribution::WeightedPolygon_t polygon;
596 polygon += std::make_pair( Vector(1.,0.,0.), 1);
597 polygon += std::make_pair( Vector(-1./3.0, 2.0*M_SQRT2/3.0,0.), 1);
598 polygon += std::make_pair( Vector(-1./3.0, -M_SQRT2/3.0, M_SQRT2/sqrt(3)), 1);
599 SphericalPointDistribution::Polygon_t newpolygon =
600 SPD.get<4>();
601 SphericalPointDistribution::Polygon_t expected = newpolygon;
602 expected.pop_front(); // remove first point
603 expected.pop_front(); // remove second point
604 expected.pop_front(); // remove third point
605 SphericalPointDistribution::Polygon_t remaining =
606 SphericalPointDistribution::matchSphericalPointDistributions(
607 polygon,
608 newpolygon);
609 CPPUNIT_ASSERT_EQUAL( expected, remaining );
610 // also slightly perturbed
611 const double amplitude = 0.05;
612 perturbPolygon(polygon, amplitude);
613 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
614 }
615
616 // test with three points, full rotation
617 {
618 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
619 SphericalPointDistribution::WeightedPolygon_t polygon;
620 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
621 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1./3.0, 2.0*M_SQRT2/3.0,0.), 47.6/180*M_PI), 1);
622 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1./3.0, -M_SQRT2/3.0, M_SQRT2/sqrt(3)), 47.6/180*M_PI), 1);
623 SphericalPointDistribution::Polygon_t newpolygon =
624 SPD.get<4>();
625 SphericalPointDistribution::Polygon_t expected = newpolygon;
626 expected.pop_front(); // remove first point
627 expected.pop_front(); // remove second point
628 expected.pop_front(); // remove third point
629 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
630 iter != expected.end(); ++iter)
631 *iter = RotationAxis.rotateVector(*iter, 47.6/180*M_PI);
632 SphericalPointDistribution::Polygon_t remaining =
633 SphericalPointDistribution::matchSphericalPointDistributions(
634 polygon,
635 newpolygon);
636 CPPUNIT_ASSERT_EQUAL( expected, remaining );
637 // also slightly perturbed
638 const double amplitude = 0.05;
639 perturbPolygon(polygon, amplitude);
640 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
641 }
642}
643
644/** UnitTest for matchSphericalPointDistributions() with four points and weights
645 * not all equal to one.
646 */
647void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_multiple()
648{
649 SphericalPointDistribution SPD(1.);
650
651 // test with four points: one point having weight of two
652 {
653 SphericalPointDistribution::WeightedPolygon_t polygon;
654 polygon += std::make_pair( Vector(1.,0.,0.), 2);
655 SphericalPointDistribution::Polygon_t newpolygon =
656 SPD.get<4>();
657 SphericalPointDistribution::Polygon_t expected;
658 expected += Vector(-0.5773502691896,-5.551115123126e-17,0.8164965809277);
659 expected += Vector(-0.5773502691896,-5.551115123126e-17,-0.8164965809277);
660 SphericalPointDistribution::Polygon_t remaining =
661 SphericalPointDistribution::matchSphericalPointDistributions(
662 polygon,
663 newpolygon);
664// std::cout << std::setprecision(13) << "Matched polygon is " << remaining << std::endl;
665// CPPUNIT_ASSERT_EQUAL( expected, remaining );
666 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, CenterAccuracy) );
667 }
668
669 // test with five points: one point having weight of two
670 {
671 SphericalPointDistribution::WeightedPolygon_t polygon;
672 polygon += std::make_pair( Vector(1.,0.,0.), 2);
673 SphericalPointDistribution::Polygon_t newpolygon =
674 SPD.get<5>();
675 SphericalPointDistribution::Polygon_t expected;
676 expected += Vector(-0.7071067811865,0.7071067811865,0);
677 expected += Vector(-0.3535533905933,-0.3535533905933,0.8660254037844);
678 expected += Vector(-0.3535533905933,-0.3535533905933,-0.8660254037844);
679 SphericalPointDistribution::Polygon_t remaining =
680 SphericalPointDistribution::matchSphericalPointDistributions(
681 polygon,
682 newpolygon);
683// std::cout << std::setprecision(13) << "Matched polygon is " << remaining << std::endl;
684// CPPUNIT_ASSERT_EQUAL( expected, remaining );
685 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, CenterAccuracy) );
686 }
687
688
689 // test with five points: one point having weight of two, one weight of one
690 {
691 SphericalPointDistribution::WeightedPolygon_t polygon;
692 polygon += std::make_pair( Vector(M_SQRT1_2,M_SQRT1_2,0.), 2);
693 polygon += std::make_pair( Vector(-1.,0.,0.), 1);
694 SphericalPointDistribution::Polygon_t newpolygon =
695 SPD.get<5>();
696 SphericalPointDistribution::Polygon_t expected;
697 expected += Vector(0.3535533786708,-0.3535533955317,-0.8660254066357);
698 expected += Vector(0.3535534025157,-0.3535533856548,0.8660254009332);
699 SphericalPointDistribution::Polygon_t remaining =
700 SphericalPointDistribution::matchSphericalPointDistributions(
701 polygon,
702 newpolygon);
703// std::cout << std::setprecision(13) << "Matched polygon is " << remaining << std::endl;
704// CPPUNIT_ASSERT_EQUAL( expected, remaining );
705 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, CenterAccuracy) );
706 }
707
708 // test with six points: two points each having weight of two
709 {
710 SphericalPointDistribution::WeightedPolygon_t polygon;
711 polygon += std::make_pair( Vector(M_SQRT1_2,-M_SQRT1_2,0.), 2);
712 polygon += std::make_pair( Vector(-M_SQRT1_2,M_SQRT1_2,0.), 2);
713 SphericalPointDistribution::Polygon_t newpolygon =
714 SPD.get<6>();
715 SphericalPointDistribution::Polygon_t expected;
716 expected += Vector(0.,0.,-1.);
717 expected += Vector(0.,0.,1.);
718 SphericalPointDistribution::Polygon_t remaining =
719 SphericalPointDistribution::matchSphericalPointDistributions(
720 polygon,
721 newpolygon);
722// std::cout << std::setprecision(13) << "Matched polygon is " << remaining << std::endl;
723// CPPUNIT_ASSERT_EQUAL( expected, remaining );
724 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, CenterAccuracy) );
725 }
726}
727
728/** UnitTest for matchSphericalPointDistributions() with five points
729 */
730void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_5()
731{
732 SphericalPointDistribution SPD(1.);
733
734 // test with one point, matching trivially
735 {
736 SphericalPointDistribution::WeightedPolygon_t polygon;
737 polygon += std::make_pair( Vector(1.,0.,0.), 1);
738 SphericalPointDistribution::Polygon_t newpolygon =
739 SPD.get<5>();
740 SphericalPointDistribution::Polygon_t expected = newpolygon;
741 expected.pop_front(); // remove first point
742 SphericalPointDistribution::Polygon_t remaining =
743 SphericalPointDistribution::matchSphericalPointDistributions(
744 polygon,
745 newpolygon);
746 CPPUNIT_ASSERT_EQUAL( expected, remaining );
747 }
748
749 // test with one point, just a flip of axis
750 {
751 SphericalPointDistribution::WeightedPolygon_t polygon;
752 polygon += std::make_pair( Vector(0.,1.,0.), 1);
753 SphericalPointDistribution::Polygon_t newpolygon =
754 SPD.get<5>();
755 SphericalPointDistribution::Polygon_t expected = newpolygon;
756 expected.pop_front(); // remove first point
757 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
758 iter != expected.end(); ++iter) {
759 std::swap((*iter)[0], (*iter)[1]);
760 (*iter)[0] *= -1.;
761 }
762 SphericalPointDistribution::Polygon_t remaining =
763 SphericalPointDistribution::matchSphericalPointDistributions(
764 polygon,
765 newpolygon);
766 CPPUNIT_ASSERT_EQUAL( expected, remaining );
767 }
768
769 // test with two points, matching trivially
770 {
771 SphericalPointDistribution::WeightedPolygon_t polygon;
772 polygon += std::make_pair( Vector(1.,0.,0.), 1);
773 polygon += std::make_pair( Vector(-1.,0.,0.), 1);
774 SphericalPointDistribution::Polygon_t newpolygon =
775 SPD.get<5>();
776 SphericalPointDistribution::Polygon_t expected = newpolygon;
777 expected.pop_front(); // remove first point
778 expected.pop_front(); // remove second point
779 SphericalPointDistribution::Polygon_t remaining =
780 SphericalPointDistribution::matchSphericalPointDistributions(
781 polygon,
782 newpolygon);
783 CPPUNIT_ASSERT_EQUAL( expected, remaining );
784 // also slightly perturbed
785 const double amplitude = 0.05;
786 perturbPolygon(polygon, amplitude);
787 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
788 }
789
790 // test with two points, full rotation
791 {
792 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
793 SphericalPointDistribution::WeightedPolygon_t polygon;
794 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180.*M_PI), 1);
795 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1.,0.,0.), 47.6/180.*M_PI), 1);
796 SphericalPointDistribution::Polygon_t newpolygon =
797 SPD.get<5>();
798 SphericalPointDistribution::Polygon_t expected = newpolygon;
799 expected.pop_front(); // remove first point
800 expected.pop_front(); // remove second point
801 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
802 iter != expected.end(); ++iter)
803 *iter = RotationAxis.rotateVector(*iter, 47.6/180.*M_PI);
804 SphericalPointDistribution::Polygon_t remaining =
805 SphericalPointDistribution::matchSphericalPointDistributions(
806 polygon,
807 newpolygon);
808 // the three remaining points sit on a plane that may be rotated arbitrarily
809 // so we cannot simply check for equality between expected and remaining
810 // hence, we just check that they are orthogonal to the first two points
811 CPPUNIT_ASSERT_EQUAL( expected.size(), remaining.size() );
812 for (SphericalPointDistribution::WeightedPolygon_t::const_iterator fixiter = polygon.begin();
813 fixiter != polygon.end(); ++fixiter) {
814 for (SphericalPointDistribution::Polygon_t::const_iterator iter = remaining.begin();
815 iter != remaining.end(); ++iter) {
816 CPPUNIT_ASSERT( (fixiter->first).IsNormalTo(*iter) );
817 }
818 }
819 }
820
821 // test with three points, matching trivially
822 {
823 SphericalPointDistribution::WeightedPolygon_t polygon;
824 polygon += std::make_pair( Vector(1.,0.,0.), 1);
825 polygon += std::make_pair( Vector(-1., 0.0, 0.0), 1);
826 polygon += std::make_pair( Vector(0.0, 1., 0.0), 1);
827 SphericalPointDistribution::Polygon_t newpolygon =
828 SPD.get<5>();
829 SphericalPointDistribution::Polygon_t expected = newpolygon;
830 expected.pop_front(); // remove first point
831 expected.pop_front(); // remove second point
832 expected.pop_front(); // remove third point
833 SphericalPointDistribution::Polygon_t remaining =
834 SphericalPointDistribution::matchSphericalPointDistributions(
835 polygon,
836 newpolygon);
837 CPPUNIT_ASSERT_EQUAL( expected, remaining );
838 // also slightly perturbed
839 const double amplitude = 0.05;
840 perturbPolygon(polygon, amplitude);
841 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
842 }
843
844 // test with three points, full rotation
845 {
846 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
847 SphericalPointDistribution::WeightedPolygon_t polygon;
848 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
849 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1., 0.0, 0.0), 47.6/180*M_PI), 1);
850 polygon += std::make_pair(RotationAxis.rotateVector(Vector(0.0, 1., 0.0), 47.6/180*M_PI), 1);
851 SphericalPointDistribution::Polygon_t newpolygon =
852 SPD.get<5>();
853 SphericalPointDistribution::Polygon_t expected = newpolygon;
854 expected.pop_front(); // remove first point
855 expected.pop_front(); // remove second point
856 expected.pop_front(); // remove third point
857 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
858 iter != expected.end(); ++iter)
859 *iter = RotationAxis.rotateVector(*iter, 47.6/180*M_PI);
860 SphericalPointDistribution::Polygon_t remaining =
861 SphericalPointDistribution::matchSphericalPointDistributions(
862 polygon,
863 newpolygon);
864 CPPUNIT_ASSERT_EQUAL( expected, remaining );
865 // also slightly perturbed
866 const double amplitude = 0.05;
867 perturbPolygon(polygon, amplitude);
868 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
869 }
870}
871
872/** UnitTest for matchSphericalPointDistributions() with six points
873 */
874void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_6()
875{
876 SphericalPointDistribution SPD(1.);
877
878 // test with one point, matching trivially
879 {
880 SphericalPointDistribution::WeightedPolygon_t polygon;
881 polygon += std::make_pair( Vector(1.,0.,0.), 1);
882 SphericalPointDistribution::Polygon_t newpolygon =
883 SPD.get<6>();
884 SphericalPointDistribution::Polygon_t expected = newpolygon;
885 expected.pop_front(); // remove first point
886 SphericalPointDistribution::Polygon_t remaining =
887 SphericalPointDistribution::matchSphericalPointDistributions(
888 polygon,
889 newpolygon);
890 CPPUNIT_ASSERT_EQUAL( expected, remaining );
891 }
892
893 // test with one point, just a flip of axis
894 {
895 SphericalPointDistribution::WeightedPolygon_t polygon;
896 polygon += std::make_pair( Vector(0.,1.,0.), 1);
897 SphericalPointDistribution::Polygon_t newpolygon =
898 SPD.get<6>();
899 SphericalPointDistribution::Polygon_t expected = newpolygon;
900 expected.pop_front(); // remove first point
901 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
902 iter != expected.end(); ++iter) {
903 std::swap((*iter)[0], (*iter)[1]);
904 (*iter)[0] *= -1.;
905 }
906 SphericalPointDistribution::Polygon_t remaining =
907 SphericalPointDistribution::matchSphericalPointDistributions(
908 polygon,
909 newpolygon);
910 CPPUNIT_ASSERT_EQUAL( expected, remaining );
911 }
912
913 // test with two points, matching trivially
914 {
915 SphericalPointDistribution::WeightedPolygon_t polygon;
916 polygon += std::make_pair( Vector(1.,0.,0.), 1);
917 polygon += std::make_pair( Vector(-1.,0.,0.), 1);
918 SphericalPointDistribution::Polygon_t newpolygon =
919 SPD.get<6>();
920 SphericalPointDistribution::Polygon_t expected = newpolygon;
921 expected.pop_front(); // remove first point
922 expected.pop_front(); // remove second spoint
923 SphericalPointDistribution::Polygon_t remaining =
924 SphericalPointDistribution::matchSphericalPointDistributions(
925 polygon,
926 newpolygon);
927 CPPUNIT_ASSERT_EQUAL( expected, remaining );
928 // also slightly perturbed
929 const double amplitude = 0.05;
930 perturbPolygon(polygon, amplitude);
931 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
932 }
933
934 // test with two points, full rotation
935 {
936 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
937 SphericalPointDistribution::WeightedPolygon_t polygon;
938 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
939 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1.,0.,0.), 47.6/180*M_PI), 1);
940 SphericalPointDistribution::Polygon_t newpolygon =
941 SPD.get<6>();
942 SphericalPointDistribution::Polygon_t expected = newpolygon;
943 expected.pop_front(); // remove first point
944 expected.pop_front(); // remove second spoint
945 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
946 iter != expected.end(); ++iter)
947 *iter = RotationAxis.rotateVector(*iter, 47.6/180*M_PI);
948 SphericalPointDistribution::Polygon_t remaining =
949 SphericalPointDistribution::matchSphericalPointDistributions(
950 polygon,
951 newpolygon);
952 // the four remaining points sit on a plane that may have been rotated arbitrarily
953 // so we cannot simply check for equality between expected and remaining
954 // hence, we just check that they are orthogonal to the first two points
955 CPPUNIT_ASSERT_EQUAL( expected.size(), remaining.size() );
956 for (SphericalPointDistribution::WeightedPolygon_t::const_iterator fixiter = polygon.begin();
957 fixiter != polygon.end(); ++fixiter) {
958 for (SphericalPointDistribution::Polygon_t::const_iterator iter = remaining.begin();
959 iter != remaining.end(); ++iter) {
960 CPPUNIT_ASSERT( (fixiter->first).IsNormalTo(*iter) );
961 }
962 }
963 }
964
965 // test with three points, matching trivially
966 {
967 SphericalPointDistribution::WeightedPolygon_t polygon;
968 polygon += std::make_pair( Vector(1.,0.,0.), 1);
969 polygon += std::make_pair( Vector(-1., 0.0, 0.0), 1);
970 polygon += std::make_pair( Vector(0.0, 1., 0.0), 1);
971 SphericalPointDistribution::Polygon_t newpolygon =
972 SPD.get<6>();
973 SphericalPointDistribution::Polygon_t expected = newpolygon;
974 expected.pop_front(); // remove first point
975 expected.pop_front(); // remove second point
976 expected.pop_front(); // remove third point
977 SphericalPointDistribution::Polygon_t remaining =
978 SphericalPointDistribution::matchSphericalPointDistributions(
979 polygon,
980 newpolygon);
981 CPPUNIT_ASSERT_EQUAL( expected, remaining );
982 // also slightly perturbed
983 const double amplitude = 0.05;
984 perturbPolygon(polygon, amplitude);
985 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
986 }
987
988 // test with three points, full rotation
989 {
990 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
991 SphericalPointDistribution::WeightedPolygon_t polygon;
992 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
993 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1., 0.0, 0.0), 47.6/180*M_PI), 1);
994 polygon += std::make_pair(RotationAxis.rotateVector(Vector(0.0, 1., 0.0), 47.6/180*M_PI), 1);
995 SphericalPointDistribution::Polygon_t newpolygon =
996 SPD.get<6>();
997 SphericalPointDistribution::Polygon_t expected = newpolygon;
998 expected.pop_front(); // remove first point
999 expected.pop_front(); // remove second point
1000 expected.pop_front(); // remove third point
1001 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
1002 iter != expected.end(); ++iter)
1003 *iter = RotationAxis.rotateVector(*iter, 47.6/180*M_PI);
1004 SphericalPointDistribution::Polygon_t remaining =
1005 SphericalPointDistribution::matchSphericalPointDistributions(
1006 polygon,
1007 newpolygon);
1008 CPPUNIT_ASSERT_EQUAL( expected, remaining );
1009 // also slightly perturbed
1010 const double amplitude = 0.05;
1011 perturbPolygon(polygon, amplitude);
1012 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
1013 }
1014}
1015
1016/** UnitTest for matchSphericalPointDistributions() with seven points
1017 */
1018void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_7()
1019{
1020 SphericalPointDistribution SPD(1.);
1021
1022 // test with one point, matching trivially
1023 {
1024 SphericalPointDistribution::WeightedPolygon_t polygon;
1025 polygon += std::make_pair( Vector(1.,0.,0.), 1);
1026 SphericalPointDistribution::Polygon_t newpolygon =
1027 SPD.get<7>();
1028 SphericalPointDistribution::Polygon_t expected = newpolygon;
1029 expected.pop_front(); // remove first point
1030 SphericalPointDistribution::Polygon_t remaining =
1031 SphericalPointDistribution::matchSphericalPointDistributions(
1032 polygon,
1033 newpolygon);
1034 CPPUNIT_ASSERT_EQUAL( expected, remaining );
1035 }
1036
1037 // test with one point, just a flip of axis
1038 {
1039 SphericalPointDistribution::WeightedPolygon_t polygon;
1040 polygon += std::make_pair( Vector(0.,1.,0.), 1);
1041 SphericalPointDistribution::Polygon_t newpolygon =
1042 SPD.get<7>();
1043 SphericalPointDistribution::Polygon_t expected = newpolygon;
1044 expected.pop_front(); // remove first point
1045 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
1046 iter != expected.end(); ++iter) {
1047 std::swap((*iter)[0], (*iter)[1]);
1048 (*iter)[0] *= -1.;
1049 }
1050 SphericalPointDistribution::Polygon_t remaining =
1051 SphericalPointDistribution::matchSphericalPointDistributions(
1052 polygon,
1053 newpolygon);
1054 CPPUNIT_ASSERT_EQUAL( expected, remaining );
1055 }
1056
1057 // test with two points, matching trivially
1058 {
1059 SphericalPointDistribution::WeightedPolygon_t polygon;
1060 polygon += std::make_pair( Vector(1.,0.,0.), 1);
1061 polygon += std::make_pair( Vector(-1.,0.,0.), 1);
1062 SphericalPointDistribution::Polygon_t newpolygon =
1063 SPD.get<7>();
1064 SphericalPointDistribution::Polygon_t expected = newpolygon;
1065 expected.pop_front(); // remove first point
1066 expected.pop_front(); // remove second point
1067 SphericalPointDistribution::Polygon_t remaining =
1068 SphericalPointDistribution::matchSphericalPointDistributions(
1069 polygon,
1070 newpolygon);
1071 CPPUNIT_ASSERT_EQUAL( expected, remaining );
1072 // also slightly perturbed
1073 const double amplitude = 0.05;
1074 perturbPolygon(polygon, amplitude);
1075 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
1076 }
1077
1078 // test with two points, full rotation
1079 {
1080 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
1081 SphericalPointDistribution::WeightedPolygon_t polygon;
1082 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
1083 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1.,0.,0.), 47.6/180*M_PI), 1);
1084 SphericalPointDistribution::Polygon_t newpolygon =
1085 SPD.get<7>();
1086 SphericalPointDistribution::Polygon_t expected = newpolygon;
1087 expected.pop_front(); // remove first point
1088 expected.pop_front(); // remove second point
1089 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
1090 iter != expected.end(); ++iter)
1091 *iter = RotationAxis.rotateVector(*iter, 47.6/180*M_PI);
1092 SphericalPointDistribution::Polygon_t remaining =
1093 SphericalPointDistribution::matchSphericalPointDistributions(
1094 polygon,
1095 newpolygon);
1096 // the five remaining points sit on a plane that may have been rotated arbitrarily
1097 // so we cannot simply check for equality between expected and remaining
1098 // hence, we just check that they are orthogonal to the first two points
1099 CPPUNIT_ASSERT_EQUAL( expected.size(), remaining.size() );
1100 for (SphericalPointDistribution::WeightedPolygon_t::const_iterator fixiter = polygon.begin();
1101 fixiter != polygon.end(); ++fixiter) {
1102 for (SphericalPointDistribution::Polygon_t::const_iterator iter = remaining.begin();
1103 iter != remaining.end(); ++iter) {
1104 CPPUNIT_ASSERT( (fixiter->first).IsNormalTo(*iter) );
1105 }
1106 }
1107 }
1108
1109 // test with three points, matching trivially
1110 {
1111 SphericalPointDistribution::WeightedPolygon_t polygon;
1112 polygon += std::make_pair( Vector(1.,0.,0.), 1);
1113 polygon += std::make_pair( Vector(-1., 0.0, 0.0), 1);
1114 polygon += std::make_pair( Vector(0.0, 1., 0.0), 1);
1115 SphericalPointDistribution::Polygon_t newpolygon =
1116 SPD.get<7>();
1117 SphericalPointDistribution::Polygon_t expected = newpolygon;
1118 expected.pop_front(); // remove first point
1119 expected.pop_front(); // remove second point
1120 expected.pop_front(); // remove third point
1121 SphericalPointDistribution::Polygon_t remaining =
1122 SphericalPointDistribution::matchSphericalPointDistributions(
1123 polygon,
1124 newpolygon);
1125 CPPUNIT_ASSERT_EQUAL( expected, remaining );
1126 // also slightly perturbed
1127 const double amplitude = 0.05;
1128 perturbPolygon(polygon, amplitude);
1129 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
1130 }
1131
1132 // test with three points, full rotation
1133 {
1134 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
1135 SphericalPointDistribution::WeightedPolygon_t polygon;
1136 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
1137 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1., 0.0, 0.0), 47.6/180*M_PI), 1);
1138 polygon += std::make_pair(RotationAxis.rotateVector(Vector(0.0, 1., 0.0), 47.6/180*M_PI), 1);
1139 SphericalPointDistribution::Polygon_t newpolygon =
1140 SPD.get<7>();
1141 SphericalPointDistribution::Polygon_t expected = newpolygon;
1142 expected.pop_front(); // remove first point
1143 expected.pop_front(); // remove second point
1144 expected.pop_front(); // remove third point
1145 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
1146 iter != expected.end(); ++iter)
1147 *iter = RotationAxis.rotateVector(*iter, 47.6/180*M_PI);
1148 SphericalPointDistribution::Polygon_t remaining =
1149 SphericalPointDistribution::matchSphericalPointDistributions(
1150 polygon,
1151 newpolygon);
1152 CPPUNIT_ASSERT_EQUAL( expected, remaining );
1153 // also slightly perturbed
1154 const double amplitude = 0.05;
1155 perturbPolygon(polygon, amplitude);
1156 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
1157 }
1158}
1159
1160/** UnitTest for matchSphericalPointDistributions() with eight points
1161 */
1162void SphericalPointDistributionTest::matchSphericalPointDistributionsTest_8()
1163{
1164 SphericalPointDistribution SPD(1.);
1165
1166 // test with one point, matching trivially
1167 {
1168 SphericalPointDistribution::WeightedPolygon_t polygon;
1169 polygon += std::make_pair( Vector(1.,0.,0.), 1);
1170 SphericalPointDistribution::Polygon_t newpolygon =
1171 SPD.get<8>();
1172 SphericalPointDistribution::Polygon_t expected = newpolygon;
1173 expected.pop_front(); // remove first point
1174 SphericalPointDistribution::Polygon_t remaining =
1175 SphericalPointDistribution::matchSphericalPointDistributions(
1176 polygon,
1177 newpolygon);
1178 CPPUNIT_ASSERT_EQUAL( expected, remaining );
1179 }
1180
1181 // test with one point, just a flip of axis
1182 {
1183 SphericalPointDistribution::WeightedPolygon_t polygon;
1184 polygon += std::make_pair( Vector(0.,1.,0.), 1);
1185 SphericalPointDistribution::Polygon_t newpolygon =
1186 SPD.get<8>();
1187 SphericalPointDistribution::Polygon_t expected = newpolygon;
1188 expected.pop_front(); // remove first point
1189 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
1190 iter != expected.end(); ++iter) {
1191 std::swap((*iter)[0], (*iter)[1]);
1192 (*iter)[0] *= -1.;
1193 }
1194 SphericalPointDistribution::Polygon_t remaining =
1195 SphericalPointDistribution::matchSphericalPointDistributions(
1196 polygon,
1197 newpolygon);
1198 CPPUNIT_ASSERT_EQUAL( expected, remaining );
1199 }
1200
1201 // test with two points, matching trivially
1202 {
1203 SphericalPointDistribution::WeightedPolygon_t polygon;
1204 polygon += std::make_pair( Vector(1.,0.,0.), 1);
1205 polygon += std::make_pair( Vector(-1.,0.,0.), 1);
1206 SphericalPointDistribution::Polygon_t newpolygon =
1207 SPD.get<8>();
1208 SphericalPointDistribution::Polygon_t expected = newpolygon;
1209 expected.pop_front(); // remove first point
1210 expected.pop_front(); // remove second point
1211 SphericalPointDistribution::Polygon_t remaining =
1212 SphericalPointDistribution::matchSphericalPointDistributions(
1213 polygon,
1214 newpolygon);
1215 CPPUNIT_ASSERT_EQUAL( expected, remaining );
1216 // also slightly perturbed
1217 const double amplitude = 0.05;
1218 perturbPolygon(polygon, amplitude);
1219 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
1220 }
1221
1222 // test with two points, full rotation
1223 {
1224 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
1225 SphericalPointDistribution::WeightedPolygon_t polygon;
1226 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
1227 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1.,0.,0.), 47.6/180*M_PI), 1);
1228 SphericalPointDistribution::Polygon_t newpolygon =
1229 SPD.get<8>();
1230 SphericalPointDistribution::Polygon_t expected = newpolygon;
1231 expected.pop_front(); // remove first point
1232 expected.pop_front(); // remove second point
1233 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
1234 iter != expected.end(); ++iter)
1235 *iter = RotationAxis.rotateVector(*iter, 47.6/180*M_PI);
1236 SphericalPointDistribution::Polygon_t remaining =
1237 SphericalPointDistribution::matchSphericalPointDistributions(
1238 polygon,
1239 newpolygon);
1240 // the six remaining points sit on two planes that may have been rotated arbitrarily
1241 // so we cannot simply check for equality between expected and remaining
1242 // hence, we just check that they are orthogonal to the first two points
1243 CPPUNIT_ASSERT_EQUAL( expected.size(), remaining.size() );
1244 for (SphericalPointDistribution::WeightedPolygon_t::const_iterator fixiter = polygon.begin();
1245 fixiter != polygon.end(); ++fixiter) {
1246 SphericalPointDistribution::Polygon_t::const_iterator expectiter = expected.begin();
1247 SphericalPointDistribution::Polygon_t::const_iterator remainiter = remaining.begin();
1248 for (;remainiter != remaining.end(); ++expectiter, ++remainiter) {
1249 // check that points in expected/remaining have same angle to the given ones
1250// CPPUNIT_ASSERT_EQUAL( (*expectiter).Angle(*fixiter), (*remainiter).Angle(*fixiter) );
1251 CPPUNIT_ASSERT( fabs( (*expectiter).Angle(fixiter->first) - (*remainiter).Angle(fixiter->first) )
1252 < std::numeric_limits<double>::epsilon()*1e4 );
1253 }
1254 }
1255 }
1256
1257 // test with three points, matching trivially
1258 {
1259 SphericalPointDistribution::WeightedPolygon_t polygon;
1260 polygon += std::make_pair( Vector(1.,0.,0.), 1);
1261 polygon += std::make_pair( Vector(-1., 0.0, 0.0), 1);
1262 polygon += std::make_pair( Vector(-1./3.0, 2.0*M_SQRT2/3.0, 0.0), 1);
1263 SphericalPointDistribution::Polygon_t newpolygon =
1264 SPD.get<8>();
1265 SphericalPointDistribution::Polygon_t expected = newpolygon;
1266 expected.pop_front(); // remove first point
1267 expected.pop_front(); // remove second point
1268 expected.pop_front(); // remove third point
1269 SphericalPointDistribution::Polygon_t remaining =
1270 SphericalPointDistribution::matchSphericalPointDistributions(
1271 polygon,
1272 newpolygon);
1273 CPPUNIT_ASSERT_EQUAL( expected, remaining );
1274 // also slightly perturbed
1275 const double amplitude = 0.05;
1276 perturbPolygon(polygon, amplitude);
1277 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
1278 }
1279
1280 // test with three points, full rotation
1281 {
1282 Line RotationAxis(zeroVec, Vector(0.2, 0.43, 0.6893248));
1283 SphericalPointDistribution::WeightedPolygon_t polygon;
1284 polygon += std::make_pair(RotationAxis.rotateVector(Vector(1.,0.,0.), 47.6/180*M_PI), 1);
1285 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1., 0.0, 0.0), 47.6/180*M_PI), 1);
1286 polygon += std::make_pair(RotationAxis.rotateVector(Vector(-1./3.0, 2.0*M_SQRT2/3.0, 0.0), 47.6/180*M_PI), 1);
1287 SphericalPointDistribution::Polygon_t newpolygon =
1288 SPD.get<8>();
1289 SphericalPointDistribution::Polygon_t expected = newpolygon;
1290 expected.pop_front(); // remove first point
1291 expected.pop_front(); // remove second point
1292 expected.pop_front(); // remove third point
1293 for (SphericalPointDistribution::Polygon_t::iterator iter = expected.begin();
1294 iter != expected.end(); ++iter)
1295 *iter = RotationAxis.rotateVector(*iter, 47.6/180*M_PI);
1296 SphericalPointDistribution::Polygon_t remaining =
1297 SphericalPointDistribution::matchSphericalPointDistributions(
1298 polygon,
1299 newpolygon);
1300 CPPUNIT_ASSERT_EQUAL( expected, remaining );
1301 // also slightly perturbed
1302 const double amplitude = 0.05;
1303 perturbPolygon(polygon, amplitude);
1304 CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
1305 }
1306}
Note: See TracBrowser for help on using the repository browser.