Ignore:
Timestamp:
Sep 12, 2016, 11:48:34 PM (9 years ago)
Author:
Frederik Heber <heber@…>
Branches:
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
Children:
2199c2
Parents:
158ecb
git-author:
Frederik Heber <heber@…> (06/12/14 07:23:12)
git-committer:
Frederik Heber <heber@…> (09/12/16 23:48:34)
Message:

Using the idea of three points giving a triangle to find rotation axis.

  • we calculate the center of either triangle and rotate the center of the ideal point distribution to match the one from the given points.
  • next we have the triangles normals as axis, take the first matching point and rotate align it.
  • we have to deal with a lot of special cases: What if only zero, one, or two points are given ...
  • in general we assume that the triangle lies relatively flat on the sphere's surface but what if the origin is in the triangle plane or even the calculated center is at the origin ...
  • TESTS: SphericalPointDistributionUnitTest working again, regression tests FragmentMolecule-cylces and StoreSaturatedFragment working.
File:
1 edited

Legend:

Unmodified
Added
Removed
  • src/Fragmentation/Exporters/unittests/SphericalPointDistributionUnitTest.cpp

    r158ecb rb67d89  
    150150}
    151151
     152void perturbPolygon(
     153    SphericalPointDistribution::Polygon_t &_polygon,
     154    double _amplitude
     155    )
     156{
     157  for (SphericalPointDistribution::Polygon_t::iterator iter = _polygon.begin();
     158      iter != _polygon.end(); ++iter) {
     159    Vector perturber;
     160    perturber.GetOneNormalVector((*iter));
     161    perturber.Scale(_amplitude);
     162    *iter = *iter + perturber;
     163    (*iter).Normalize();
     164  }
     165}
     166
     167static
     168bool areEqualToWithinBounds(
     169    const SphericalPointDistribution::Polygon_t &_polygon,
     170    const SphericalPointDistribution::Polygon_t &_otherpolygon,
     171    double _amplitude
     172    )
     173{
     174  // same size?
     175  if (_polygon.size() != _otherpolygon.size())
     176    return false;
     177  // same points ? We just check witrh trivial mapping, nothing fancy ...
     178  bool status = true;
     179  SphericalPointDistribution::Polygon_t::const_iterator iter = _polygon.begin();
     180  SphericalPointDistribution::Polygon_t::const_iterator otheriter = _otherpolygon.begin();
     181  for (; iter != _polygon.end(); ++iter, ++otheriter) {
     182    status &= (*iter - *otheriter).Norm() < _amplitude;
     183  }
     184  return status;
     185}
     186
     187/** UnitTest for areEqualToWithinBounds()
     188 */
     189void SphericalPointDistributionTest::areEqualToWithinBoundsTest()
     190{
     191  // test with no points
     192  {
     193    SphericalPointDistribution::Polygon_t polygon;
     194    SphericalPointDistribution::Polygon_t expected = polygon;
     195    CPPUNIT_ASSERT( areEqualToWithinBounds(polygon, expected, std::numeric_limits<double>::epsilon()*1e2) );
     196  }
     197  // test with one point
     198  {
     199    SphericalPointDistribution::Polygon_t polygon;
     200    polygon += Vector(1.,0.,0.);
     201    SphericalPointDistribution::Polygon_t expected = polygon;
     202    CPPUNIT_ASSERT( areEqualToWithinBounds(polygon, expected, std::numeric_limits<double>::epsilon()*1e2) );
     203  }
     204  // test with two points
     205  {
     206    SphericalPointDistribution::Polygon_t polygon;
     207    polygon += Vector(1.,0.,0.);
     208    polygon += Vector(0.,1.,0.);
     209    SphericalPointDistribution::Polygon_t expected = polygon;
     210    CPPUNIT_ASSERT( areEqualToWithinBounds(polygon, expected, std::numeric_limits<double>::epsilon()*1e2) );
     211  }
     212
     213  // test with two points in different order: THIS GOES WRONG: We only check trivially
     214  {
     215    SphericalPointDistribution::Polygon_t polygon;
     216    polygon += Vector(1.,0.,0.);
     217    polygon += Vector(0.,1.,0.);
     218    SphericalPointDistribution::Polygon_t expected;
     219    expected += Vector(0.,1.,0.);
     220    expected += Vector(1.,0.,0.);
     221    CPPUNIT_ASSERT( !areEqualToWithinBounds(polygon, expected, std::numeric_limits<double>::epsilon()*1e2) );
     222  }
     223
     224  // test with two different points
     225  {
     226    SphericalPointDistribution::Polygon_t polygon;
     227    polygon += Vector(1.,0.,0.);
     228    polygon += Vector(0.,1.,0.);
     229    SphericalPointDistribution::Polygon_t expected;
     230    expected += Vector(1.01,0.,0.);
     231    expected += Vector(0.,1.,0.);
     232    CPPUNIT_ASSERT( areEqualToWithinBounds(polygon, expected, 0.05) );
     233    CPPUNIT_ASSERT( !areEqualToWithinBounds(polygon, expected, 0.005) );
     234  }
     235
     236  // test with different number of points
     237  {
     238    SphericalPointDistribution::Polygon_t polygon;
     239    polygon += Vector(1.,0.,0.);
     240    polygon += Vector(0.,1.,0.);
     241    SphericalPointDistribution::Polygon_t expected;
     242    expected += Vector(0.,1.,0.);
     243    CPPUNIT_ASSERT( !areEqualToWithinBounds(polygon, expected, 0.05) );
     244  }
     245}
     246
     247
    152248/** UnitTest for matchSphericalPointDistributions() with three points
    153249 */
     
    205301            newpolygon);
    206302    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     303    // also slightly perturbed
     304    const double amplitude = 0.05;
     305    perturbPolygon(polygon, amplitude);
     306    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    207307  }
    208308
     
    227327            newpolygon);
    228328    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     329    // also slightly perturbed
     330    const double amplitude = 0.05;
     331    perturbPolygon(polygon, amplitude);
     332    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    229333  }
    230334
     
    241345            newpolygon);
    242346    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     347    // also slightly perturbed
     348    const double amplitude = 0.05;
     349    perturbPolygon(polygon, amplitude);
     350    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    243351  }
    244352
     
    260368            newpolygon);
    261369    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     370    // also slightly perturbed
     371    const double amplitude = 0.05;
     372    perturbPolygon(polygon, amplitude);
     373    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    262374  }
    263375}
     
    318430            newpolygon);
    319431    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     432    // also slightly perturbed
     433    const double amplitude = 0.05;
     434    perturbPolygon(polygon, amplitude);
     435    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
     436  }
     437
     438  // test with two points, matching trivially, also with slightly perturbed
     439  {
     440    SphericalPointDistribution::Polygon_t polygon;
     441    polygon += Vector(1.,0.,0.), Vector(-1./3.0, 2.0*M_SQRT2/3.0,0.);
     442    SphericalPointDistribution::Polygon_t newpolygon =
     443        SPD.get<4>();
     444    SphericalPointDistribution::Polygon_t expected = newpolygon;
     445    expected.pop_front(); // remove first point
     446    expected.pop_front(); // remove second point
     447    SphericalPointDistribution::Polygon_t remaining =
     448        SphericalPointDistribution::matchSphericalPointDistributions(
     449            polygon,
     450            newpolygon);
     451    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     452    // also slightly perturbed
     453    const double amplitude = 0.05;
     454    perturbPolygon(polygon, amplitude);
     455    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    320456  }
    321457
     
    340476            newpolygon);
    341477    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) );
    342482  }
    343483
     
    360500            newpolygon);
    361501    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     502    // also slightly perturbed
     503    const double amplitude = 0.05;
     504    perturbPolygon(polygon, amplitude);
     505    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    362506  }
    363507
     
    384528            newpolygon);
    385529    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     530    // also slightly perturbed
     531    const double amplitude = 0.05;
     532    perturbPolygon(polygon, amplitude);
     533    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    386534  }
    387535}
     
    442590            newpolygon);
    443591    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     592    // also slightly perturbed
     593    const double amplitude = 0.05;
     594    perturbPolygon(polygon, amplitude);
     595    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    444596  }
    445597
     
    494646            newpolygon);
    495647    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     648    // also slightly perturbed
     649    const double amplitude = 0.05;
     650    perturbPolygon(polygon, amplitude);
     651    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    496652  }
    497653
     
    518674            newpolygon);
    519675    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     676    // also slightly perturbed
     677    const double amplitude = 0.05;
     678    perturbPolygon(polygon, amplitude);
     679    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    520680  }
    521681}
     
    576736            newpolygon);
    577737    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     738    // also slightly perturbed
     739    const double amplitude = 0.05;
     740    perturbPolygon(polygon, amplitude);
     741    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    578742  }
    579743
     
    628792            newpolygon);
    629793    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     794    // also slightly perturbed
     795    const double amplitude = 0.05;
     796    perturbPolygon(polygon, amplitude);
     797    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    630798  }
    631799
     
    652820            newpolygon);
    653821    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     822    // also slightly perturbed
     823    const double amplitude = 0.05;
     824    perturbPolygon(polygon, amplitude);
     825    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    654826  }
    655827}
     
    710882            newpolygon);
    711883    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     884    // also slightly perturbed
     885    const double amplitude = 0.05;
     886    perturbPolygon(polygon, amplitude);
     887    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    712888  }
    713889
     
    762938            newpolygon);
    763939    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     940    // also slightly perturbed
     941    const double amplitude = 0.05;
     942    perturbPolygon(polygon, amplitude);
     943    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    764944  }
    765945
     
    786966            newpolygon);
    787967    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     968    // also slightly perturbed
     969    const double amplitude = 0.05;
     970    perturbPolygon(polygon, amplitude);
     971    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    788972  }
    789973}
     
    8441028            newpolygon);
    8451029    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     1030    // also slightly perturbed
     1031    const double amplitude = 0.05;
     1032    perturbPolygon(polygon, amplitude);
     1033    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    8461034  }
    8471035
     
    9001088            newpolygon);
    9011089    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     1090    // also slightly perturbed
     1091    const double amplitude = 0.05;
     1092    perturbPolygon(polygon, amplitude);
     1093    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    9021094  }
    9031095
     
    9241116            newpolygon);
    9251117    CPPUNIT_ASSERT_EQUAL( expected, remaining );
     1118    // also slightly perturbed
     1119    const double amplitude = 0.05;
     1120    perturbPolygon(polygon, amplitude);
     1121    CPPUNIT_ASSERT( areEqualToWithinBounds(expected, remaining, amplitude) );
    9261122  }
    9271123}
Note: See TracChangeset for help on using the changeset viewer.