From 044e49906f4a8b87b2d950f2692cc4b3cfc59887 Mon Sep 17 00:00:00 2001
From: Michael Wedel <michael.wedel@psi.ch>
Date: Mon, 4 Aug 2014 16:55:09 +0200
Subject: [PATCH] Refs #9993. Completing work on PointGroup

Added documentation and some more tests.
---
 .../inc/MantidGeometry/Crystal/PointGroup.h   |  10 +-
 .../Geometry/src/Crystal/PointGroup.cpp       | 114 ++++++++++++++----
 .../src/Crystal/SymmetryOperation.cpp         |   2 +-
 .../Framework/Geometry/test/PointGroupTest.h  |  96 +++++++++++++--
 4 files changed, 183 insertions(+), 39 deletions(-)

diff --git a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Crystal/PointGroup.h b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Crystal/PointGroup.h
index 875189c6344..1f29aabb2b0 100644
--- a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Crystal/PointGroup.h
+++ b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Crystal/PointGroup.h
@@ -31,18 +31,22 @@ namespace Geometry
     /// Return true if the hkls are in same group
     virtual bool isEquivalent(Kernel::V3D hkl, Kernel::V3D hkl2) = 0;
 
+    /// Returns a vector with all equivalent hkls
     std::vector<Kernel::V3D> getEquivalents(const Kernel::V3D &hkl) const;
-    std::set<Kernel::V3D> getEquivalentSet(const Kernel::V3D &hkl) const;
-
+    /// Returns the same hkl for all equivalent hkls
     Kernel::V3D getReflectionFamily(const Kernel::V3D &hkl) const;
 
   protected:
     PointGroup();
 
     void addSymmetryOperation(const SymmetryOperation_const_sptr &symmetryOperation);
-    void calculateTransformationMatrices(const std::vector<SymmetryOperation_const_sptr> &symmetryOperations);
     std::vector<SymmetryOperation_const_sptr> getSymmetryOperations() const;
 
+    std::vector<Kernel::IntMatrix> generateTransformationMatrices(const std::vector<SymmetryOperation_const_sptr> &symmetryOperations);
+    void setTransformationMatrices(const std::vector<Kernel::IntMatrix> &matrices);
+
+    std::set<Kernel::V3D> getEquivalentSet(const Kernel::V3D &hkl) const;
+
     std::vector<SymmetryOperation_const_sptr> m_symmetryOperations;
     std::vector<Kernel::IntMatrix> m_transformationMatrices;
   };
diff --git a/Code/Mantid/Framework/Geometry/src/Crystal/PointGroup.cpp b/Code/Mantid/Framework/Geometry/src/Crystal/PointGroup.cpp
index c596acef568..db82f8e36f3 100644
--- a/Code/Mantid/Framework/Geometry/src/Crystal/PointGroup.cpp
+++ b/Code/Mantid/Framework/Geometry/src/Crystal/PointGroup.cpp
@@ -11,6 +11,20 @@ namespace Geometry
   using Kernel::V3D;
   using Kernel::IntMatrix;
 
+  /**
+   * Returns all equivalent reflections for the supplied hkl.
+   *
+   * This method returns a vector containing all equivalent hkls for the supplied one.
+   * It depends on the internal state of the pointgroup object (e.g. which symmetry operations
+   * and therefore, which transformation matrices are present). This internal state
+   * is unique for each concrete point group and is set in the constructor.
+   *
+   * The returned vector always contains a set of unique hkls, so for special hkls like (100),
+   * it has fewer entries than for a general hkl. See also PointGroup::getEquivalentSet.
+   *
+   * @param hkl :: Arbitrary hkl
+   * @return :: std::vector containing all equivalent hkls.
+   */
   std::vector<V3D> PointGroup::getEquivalents(const V3D &hkl) const
   {
     std::set<V3D> equivalents = getEquivalentSet(hkl);
@@ -18,17 +32,45 @@ namespace Geometry
     return std::vector<V3D>(equivalents.rbegin(), equivalents.rend());
   }
 
+  /**
+   * Returns the same V3D for all equivalent hkls.
+   *
+   * This method is closely related to PointGroup::getEquivalents. It returns the same
+   * V3D for all hkls of one "family". For example in a cubic point group it will return (100)
+   * for (001), (010), (0-10), etc.
+   *
+   * It can be used to generate a set of symmetry independent hkls, useful for example
+   * in powder diffraction.
+   *
+   * @param hkl :: Arbitrary hkl
+   * @return :: hkl specific to a family of index-triplets
+   */
   V3D PointGroup::getReflectionFamily(const Kernel::V3D &hkl) const
   {
     return *getEquivalentSet(hkl).rbegin();
   }
 
+  /// Protected constructor - can not be used directly.
   PointGroup::PointGroup() :
       m_symmetryOperations(),
       m_transformationMatrices()
   {
   }
 
+  /**
+   * Generates a set of hkls
+   *
+   * This method applies all transformation matrices to the supplied hkl and puts
+   * it into a set, which is returned in the end. Using a set ensures that each hkl
+   * occurs once and only once. This set is the set of equivalent hkls, specific to
+   * a concrete point group.
+   *
+   * The symmetry operations need to be set prior to calling this method by a call
+   * to PointGroup::setTransformationMatrices.
+   *
+   * @param hkl :: Arbitrary hkl
+   * @return :: set of hkls.
+   */
   std::set<V3D> PointGroup::getEquivalentSet(const Kernel::V3D &hkl) const
   {
       std::set<V3D> equivalents;
@@ -42,22 +84,48 @@ namespace Geometry
       return equivalents;
   }
 
+  /// Adds a symmetry operation to the point group.
   void PointGroup::addSymmetryOperation(const SymmetryOperation_const_sptr &symmetryOperation)
   {
       m_symmetryOperations.push_back(symmetryOperation);
   }
 
-  void PointGroup::calculateTransformationMatrices(const std::vector<SymmetryOperation_const_sptr> &symmetryOperations)
+  /// Returns all symmetry operations stored in the point group.
+  std::vector<SymmetryOperation_const_sptr> PointGroup::getSymmetryOperations() const
   {
-      m_transformationMatrices.clear();
+      return m_symmetryOperations;
+  }
 
-      std::vector<IntMatrix> trans;
-      trans.push_back(IntMatrix(3,3,true));
+  /**
+   * Returns all transformation matrices generated by a list of symmetry operations
+   *
+   * This method takes a vector of symmetry operations and returns the resulting set of
+   * transformation matrices. It does so by applying the first symmetry operation (order - 1) times
+   * to an identity matrix which results in a list of (order - 1) matrices. Then it applies
+   * the second operation to all these matrices, and so on.
+   *
+   * Using this method, all point groups can be described using a maximum of four
+   * symmetry operations. m-3m for example, which is defined in PointGroupLaue13,
+   * needs only four symmetry operations to generate all 48 transformation matrices.
+   * It does not matter which symmetry operations are chosen, as long
+   * as the chosen set generates all (but not more than) present in the point group. For 2/m, one
+   * could for example either choose (m and 2) or (m and -1) or (2 and -1).
+   *
+   * All 32 point groups can be found in:
+   *    International Tables for Crystallography A, 2006, pp. 770 - 790 (Table 10.1.2.2)
+   *
+   * @param symmetryOperations
+   * @return
+   */
+  std::vector<Kernel::IntMatrix> PointGroup::generateTransformationMatrices(const std::vector<SymmetryOperation_const_sptr> &symmetryOperations)
+  {
+      std::vector<IntMatrix> matrices;
+      matrices.push_back(IntMatrix(3,3,true));
 
       for(std::vector<SymmetryOperation_const_sptr>::const_iterator symOp = symmetryOperations.begin();
           symOp != symmetryOperations.end();
           ++symOp) {
-          std::vector<IntMatrix> currentMatrices(trans);
+          std::vector<IntMatrix> currentMatrices(matrices);
 
           for(std::vector<IntMatrix>::const_iterator currentMatrix = currentMatrices.begin();
               currentMatrix != currentMatrices.end();
@@ -65,24 +133,26 @@ namespace Geometry
               IntMatrix transformed = *currentMatrix;
               for(size_t i = 0; i < (*symOp)->order() - 1; ++i) {
                   transformed = (*symOp)->apply(transformed);
-                  trans.push_back(transformed);
+                  matrices.push_back(transformed);
               }
           }
       }
 
-      m_transformationMatrices = std::vector<IntMatrix>(trans.begin(), trans.end());
+      return matrices;
   }
 
-  std::vector<SymmetryOperation_const_sptr> PointGroup::getSymmetryOperations() const
+  /// Sets the transformation matrices.
+  void PointGroup::setTransformationMatrices(const std::vector<Kernel::IntMatrix> &matrices)
   {
-      return m_symmetryOperations;
+      m_transformationMatrices = matrices;
   }
 
+
   PointGroupLaue1::PointGroupLaue1()
   {
       addSymmetryOperation(boost::make_shared<const SymOpInversion>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue1::getName()
@@ -104,7 +174,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpRotationTwoFoldY>());
       addSymmetryOperation(boost::make_shared<const SymOpMirrorPlaneY>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue2::getName()
@@ -126,7 +196,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpRotationTwoFoldZ>());
       addSymmetryOperation(boost::make_shared<const SymOpMirrorPlaneZ>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue3::getName()
@@ -149,7 +219,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpRotationTwoFoldY>());
       addSymmetryOperation(boost::make_shared<const SymOpMirrorPlaneZ>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue4::getName()
@@ -173,7 +243,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpRotationFourFoldZ>());
       addSymmetryOperation(boost::make_shared<const SymOpMirrorPlaneZ>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue5::getName()
@@ -198,7 +268,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpRotationTwoFoldX>());
       addSymmetryOperation(boost::make_shared<const SymOpMirrorPlaneZ>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue6::getName()
@@ -225,7 +295,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpRotationThreeFoldZHexagonal>());
       addSymmetryOperation(boost::make_shared<const SymOpInversion>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue7::getName()
@@ -249,7 +319,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpInversion>());
       addSymmetryOperation(boost::make_shared<const SymOpMirrorPlane210Hexagonal>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue8::getName()
@@ -275,7 +345,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpInversion>());
       addSymmetryOperation(boost::make_shared<const SymOpRotationTwoFold210Hexagonal>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue9::getName()
@@ -300,7 +370,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpRotationSixFoldZHexagonal>());
       addSymmetryOperation(boost::make_shared<const SymOpInversion>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue10::getName()
@@ -326,7 +396,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpRotationTwoFoldXHexagonal>());
       addSymmetryOperation(boost::make_shared<const SymOpMirrorPlaneZ>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue11::getName()
@@ -357,7 +427,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpMirrorPlaneY>());
       addSymmetryOperation(boost::make_shared<const SymOpInversion>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue12::getName()
@@ -388,7 +458,7 @@ namespace Geometry
       addSymmetryOperation(boost::make_shared<const SymOpMirrorPlaneY>());
       addSymmetryOperation(boost::make_shared<const SymOpInversion>());
 
-      calculateTransformationMatrices(getSymmetryOperations());
+      setTransformationMatrices(generateTransformationMatrices(getSymmetryOperations()));
   }
 
   std::string PointGroupLaue13::getName()
diff --git a/Code/Mantid/Framework/Geometry/src/Crystal/SymmetryOperation.cpp b/Code/Mantid/Framework/Geometry/src/Crystal/SymmetryOperation.cpp
index dcfed778a5b..820c8757f7e 100644
--- a/Code/Mantid/Framework/Geometry/src/Crystal/SymmetryOperation.cpp
+++ b/Code/Mantid/Framework/Geometry/src/Crystal/SymmetryOperation.cpp
@@ -170,7 +170,7 @@ SymOpRotationThreeFold111::SymOpRotationThreeFold111() :
     setMatrixFromArray(rotThreeFold111);
 }
 
-/* 3-fold rotation axes */
+/* 6-fold rotation axes */
 /// 6-fold rotation around z-axis, hexagonal coordinate system
 SymOpRotationSixFoldZHexagonal::SymOpRotationSixFoldZHexagonal() :
     SymmetryOperation(6, Kernel::IntMatrix(3, 3), "6 [001]h")
diff --git a/Code/Mantid/Framework/Geometry/test/PointGroupTest.h b/Code/Mantid/Framework/Geometry/test/PointGroupTest.h
index d88793cedfe..6967b29ec58 100644
--- a/Code/Mantid/Framework/Geometry/test/PointGroupTest.h
+++ b/Code/Mantid/Framework/Geometry/test/PointGroupTest.h
@@ -2,6 +2,8 @@
 #define MANTID_GEOMETRY_POINTGROUPTEST_H_
 
 #include <cxxtest/TestSuite.h>
+#include <gmock/gmock.h>
+
 #include "MantidKernel/Timer.h"
 #include "MantidKernel/Strings.h"
 #include "MantidKernel/System.h"
@@ -25,8 +27,10 @@ public:
       if (pgs[i]->getName().substr(0, name.size()) == name)
       {
         std::vector<V3D> equivalents = pgs[i]->getEquivalents(hkl);
+        // check that the number of equivalent reflections is as expected.
+        TSM_ASSERT_EQUALS(name + ": Expected " + std::to_string(numEquiv) + " equivalents, got " + std::to_string(equivalents.size()) + " instead.", equivalents.size(), numEquiv);
 
-        TS_ASSERT_EQUALS(equivalents.size(), numEquiv);
+        // get reflection family for this hkl
         V3D family = pgs[i]->getReflectionFamily(hkl);
 
         for (size_t j=0; j<numEquiv; j++)
@@ -36,7 +40,10 @@ public:
           {
             TSM_ASSERT( name + " : " + hkl.toString() + " is not equivalent to " +  equiv[j].toString(), false);
           }
+
+          // make sure family for equiv[j] is the same as the one for hkl
           TS_ASSERT_EQUALS(pgs[i]->getReflectionFamily(equiv[j]), family);
+          // also make sure that current equivalent is in the collection of equivalents.
           TS_ASSERT_DIFFERS(std::find(equivalents.begin(), equivalents.end(), equiv[j]), equivalents.end());
         }
 
@@ -108,24 +115,87 @@ public:
     }
   }
 
-  void testMany()
+  void testConstruction()
+  {
+      TestablePointGroup defaultPointgroup;
+
+      TS_ASSERT_EQUALS(defaultPointgroup.m_symmetryOperations.size(), 0);
+      TS_ASSERT_EQUALS(defaultPointgroup.m_transformationMatrices.size(), 0);
+  }
+
+  void testAddSymmetryOperation()
   {
-      PointGroupLaue13 m3m;
-      V3D testHKL(1, -1, 1);
+      TestablePointGroup pg;
 
-      for(size_t i = 0; i < 1e6; ++i) {
-          //std::vector<V3D> e = m3m.getEquivalents(testHKL);
-          //V3D t = *e.begin();
+      TS_ASSERT_EQUALS(pg.getSymmetryOperations().size(), 0);
 
-          V3D family = m3m.getReflectionFamily(testHKL);
-          double x = family.X() + 1.0;
+      SymmetryOperation_const_sptr symOp(new SymOpInversion);
+      pg.addSymmetryOperation(symOp);
 
-//          for(size_t j = 0; j < e.size(); ++j) {
-//              std::cout << e[j] << std::endl;
-//          }
-      }
+      std::vector<SymmetryOperation_const_sptr> ops = pg.getSymmetryOperations();
+
+      TS_ASSERT_EQUALS(ops.size(), 1);
+      TS_ASSERT_EQUALS(ops[0], symOp);
   }
 
+  void testSetTransformationMatrices()
+  {
+      TestablePointGroup pg;
+
+      std::vector<IntMatrix> matrices(1, IntMatrix(3, 3, true));
+      pg.setTransformationMatrices(matrices);
+
+      TS_ASSERT_EQUALS(pg.m_transformationMatrices.size(), 1);
+      TS_ASSERT_EQUALS(pg.m_transformationMatrices[0], IntMatrix(3, 3, true));
+  }
+
+  void testGenerateTransformationMatrices()
+  {
+      TestablePointGroup pg;
+
+      SymmetryOperation_const_sptr identity(new SymOpIdentity);
+      SymmetryOperation_const_sptr inversion(new SymOpInversion);
+      SymmetryOperation_const_sptr mirror(new SymOpMirrorPlaneZ);
+      SymmetryOperation_const_sptr twoFold(new SymOpRotationTwoFoldZ);
+
+      pg.addSymmetryOperation(mirror);
+      pg.addSymmetryOperation(twoFold);
+
+      std::vector<SymmetryOperation_const_sptr> ops = pg.getSymmetryOperations();
+      TS_ASSERT_EQUALS(ops.size(), 2);
+
+      std::vector<IntMatrix> matrices = pg.generateTransformationMatrices(ops);
+
+      // Mirror and 2-fold axis generate inversion, identity is implicit.
+      TS_ASSERT_EQUALS(matrices.size(), 4);
+
+      auto matrixVectorBegin = matrices.begin();
+      auto matrixVectorEnd = matrices.end();
+
+      IntMatrix identityMatrix(3, 3, true);
+
+      TS_ASSERT_DIFFERS(std::find(matrixVectorBegin, matrixVectorEnd, identity->apply(identityMatrix)), matrixVectorEnd);
+      TS_ASSERT_DIFFERS(std::find(matrixVectorBegin, matrixVectorEnd, inversion->apply(identityMatrix)), matrixVectorEnd);
+      TS_ASSERT_DIFFERS(std::find(matrixVectorBegin, matrixVectorEnd, mirror->apply(identityMatrix)), matrixVectorEnd);
+      TS_ASSERT_DIFFERS(std::find(matrixVectorBegin, matrixVectorEnd, twoFold->apply(identityMatrix)), matrixVectorEnd);
+
+      TS_ASSERT_DIFFERS(matrices[0], matrices[1]);
+  }
+
+private:
+  class TestablePointGroup : public PointGroup
+  {
+      friend class PointGroupTest;
+
+  public:
+      TestablePointGroup() : PointGroup()
+      { }
+      ~TestablePointGroup() {}
+
+      MOCK_METHOD0(getName, std::string());
+      MOCK_METHOD2(isEquivalent, bool(V3D hkl, V3D hkl2));
+  };
+
 };
 
 
-- 
GitLab