Skip to content
Snippets Groups Projects
Commit 216c082d authored by Michael Wedel's avatar Michael Wedel
Browse files

Refs #13537. Finalizing HKLGenerator

parent 6e313991
No related branches found
No related tags found
No related merge requests found
......@@ -10,7 +10,69 @@
namespace Mantid {
namespace Geometry {
/** HKLGenerator : TODO: DESCRIPTION
/** HKLGenerator
HKLGenerator is a pseudo-container that helps in generating actual
containers with V3D-objects representing Miller indices (HKL).
It's a common task to generate lists of Miller indices. The simplest
way of doing that is to simply create a nested loop structure that
goes through all combinations of H, K, L within some limits and then
put them into a container (vector, list, set, ...):
for(int h = -hmin; h <= hmax; ++h) {
for(int k = -kmin; k <= kmax; ++k) {
for(int l = -lmin; l <= lmax; ++l) {
hkls.push_back(V3D(h, k, l));
}
}
}
In most cases the list is not used like that, instead it's filtered
using some criteria, for example a certain range of d-spacings or
others, so the reflection needs to be checked first:
...
hkl = V3D(h, k, l)
if(isOk(hkl)) {
hkls.push_back(hkl);
}
...
Instead of explicitly stating the triple-loop, HKLGenerator provides
a shorter way for this process using a const_iterator. The first code
example then becomes this:
HKLGenerator generator(V3D(hmin, kmin, lmin), V3D(hmax, kmax, lmax));
for(auto hkl = generator.begin(); hkl != generator.end(); ++hkl) {
hkls.push_back(*hkl);
}
Or even shorter, using std::copy:
HKLGenerator generator(V3D(hmin, kmin, lmin), V3D(hmax, kmax, lmax));
std::copy(generator.begin(), generator.end(), std::back_inserter(hkls));
It's also possible to use filters this way, but in a much easier fashion,
using std::copy_remove_if (pre C++11) or std::copy_if (C++11):
// pre C++11
std::copy_remove_if(generator.begin(), generator.end(),
std::back_inserter(hkls), isNotOk)
// C++11
std::copy_if(generator.begin(), generator.end(), std::back_inserter(hkls),
isOk)
See the documentation for HKLFilter for more details on how to perform
actual filtering.
Please be aware that the iterator increments infinitely if it passes the
specified maximimum HKL. In that case K and L remain constant while H is
incremented (until it overflows).
@author Michael Wedel, ESS
@date 23/09/2015
Copyright &copy; 2015 ISIS Rutherford Appleton Laboratory, NScD Oak Ridge
National Laboratory & European Spallation Source
......@@ -35,10 +97,23 @@ namespace Geometry {
*/
class MANTID_GEOMETRY_DLL HKLGenerator {
public:
/**
* @brief The const_iterator class
*
* This class uses boost::iterator_facade to provide a const_iterator
* interface to HKLGenerator. This makes it very easy to use the standard
* library algorithms. It's a forward iterator, so it's not possible to use
* operator[] for random access of specific HKLs (this would have to be
* done on the generated container if that supports it).
*
* While the iterators can be instantiated directly, the intention is to
* use HKLGenerator::begin() and HKLGenerator::end().
*/
class const_iterator
: public boost::iterator_facade<const_iterator, const Kernel::V3D &,
boost::forward_traversal_tag> {
public:
/// Default constructor, requirement from boost::iterator_facade
const_iterator() {}
explicit const_iterator(const Kernel::V3D &current);
......@@ -46,21 +121,20 @@ public:
explicit const_iterator(const Kernel::V3D &hklMin,
const Kernel::V3D &hklMax);
explicit const_iterator(const Kernel::V3D &hklMin,
const Kernel::V3D &hklMax,
const Kernel::V3D &current);
private:
// Required for boost::iterator_facade to work
friend class boost::iterator_core_access;
void increment();
bool equal(const const_iterator &other) const {
/// Returns true if other is at the same position
inline bool equal(const const_iterator &other) const {
return this->m_h == other.m_h && this->m_k == other.m_k &&
this->m_l == other.m_l;
}
const Kernel::V3D &dereference() const { return m_hkl; }
/// Returns a const reference to the currently pointed at HKL.
inline const Kernel::V3D &dereference() const { return m_hkl; }
int m_h, m_k, m_l;
Kernel::V3D m_hkl;
......@@ -77,19 +151,29 @@ public:
virtual ~HKLGenerator() {}
size_t size() const { return m_size; }
/// Returns the number of HKLs to be generated.
inline size_t size() const { return m_size; }
const_iterator begin() const;
const_iterator end() const;
/// Returns an iterator to the beginning of the sequence.
inline const const_iterator &begin() const { return m_begin; }
/// Returns an iterator which "points at" one element past the end.
inline const const_iterator &end() const { return m_end; }
private:
Kernel::V3D getEndHKL() const;
size_t getSize(const Kernel::V3D &min, const Kernel::V3D &max) const;
const_iterator getBeginIterator() const;
const_iterator getEndIterator() const;
Kernel::V3D getEndHKL() const;
Kernel::V3D m_hklMin;
Kernel::V3D m_hklMax;
size_t m_size;
const_iterator m_begin;
const_iterator m_end;
};
} // namespace Geometry
......
......@@ -5,79 +5,91 @@ namespace Geometry {
using namespace Kernel;
/// Constructs a generator that creates all indices from hklMin to hklMax.
HKLGenerator::HKLGenerator(const Kernel::V3D &hklMin, const Kernel::V3D &hklMax)
: m_hklMin(hklMin), m_hklMax(hklMax), m_size(getSize(m_hklMin, m_hklMax)) {}
: m_hklMin(hklMin), m_hklMax(hklMax), m_size(getSize(m_hklMin, m_hklMax)),
m_begin(getBeginIterator()), m_end(getEndIterator()) {}
/// Constructs a generator that creates all indices from -hklMinMax to
/// hklMinMax.
HKLGenerator::HKLGenerator(const Kernel::V3D &hklMinMax)
: m_hklMin(hklMinMax * -1), m_hklMax(hklMinMax),
m_size(getSize(m_hklMin, m_hklMax)) {}
m_size(getSize(m_hklMin, m_hklMax)), m_begin(getBeginIterator()),
m_end(getEndIterator()) {}
/// Constructs a generator that creates all indices from -h,-k,-l to h,k,l.
HKLGenerator::HKLGenerator(int hMinMax, int kMinMax, int lMinMax) {
m_hklMax = V3D(hMinMax, kMinMax, lMinMax);
m_hklMin = m_hklMax * -1;
m_size = getSize(m_hklMin, m_hklMax);
m_begin = getBeginIterator();
m_end = getEndIterator();
}
/// Constructs a generator that creates all indices for the given cell up to
/// dMin.
HKLGenerator::HKLGenerator(const UnitCell &unitCell, double dMin) {
m_hklMax = V3D(unitCell.a() / dMin, unitCell.b() / dMin, unitCell.c() / dMin);
m_hklMin = m_hklMax * -1;
m_size = getSize(m_hklMin, m_hklMax);
}
HKLGenerator::const_iterator HKLGenerator::end() const {
return HKLGenerator::const_iterator(getEndHKL());
}
V3D HKLGenerator::getEndHKL() const {
return V3D(m_hklMax.X() + 1, m_hklMin.Y(), m_hklMin.Z());
m_begin = getBeginIterator();
m_end = getEndIterator();
}
/// Returns the number of indices between min and max.
size_t HKLGenerator::getSize(const V3D &min, const V3D &max) const {
V3D diff = (max - min) + V3D(1, 1, 1);
return static_cast<size_t>(diff.X() * diff.Y() * diff.Z());
}
HKLGenerator::const_iterator HKLGenerator::begin() const {
/// Constructs an iterator that points to the beginning of the sequence.
HKLGenerator::const_iterator HKLGenerator::getBeginIterator() const {
return HKLGenerator::const_iterator(m_hklMin, m_hklMax);
}
/// Constructs an iterator that points to an HKL one past the maximum.
HKLGenerator::const_iterator HKLGenerator::getEndIterator() const {
return HKLGenerator::const_iterator(getEndHKL());
}
/// Returns the HKL "one past the maximum".
V3D HKLGenerator::getEndHKL() const {
return V3D(m_hklMax.X() + 1, m_hklMin.Y(), m_hklMin.Z());
}
/// Return an iterator with min = max = current.
HKLGenerator::const_iterator::const_iterator(const V3D &current)
: m_h(static_cast<int>(current.X())), m_k(static_cast<int>(current.Y())),
m_l(static_cast<int>(current.Z())), m_hMin(m_h), m_hMax(m_h), m_kMin(m_k),
m_kMax(m_k), m_lMin(m_l), m_lMax(m_l) {}
m_l(static_cast<int>(current.Z())), m_hkl(current), m_hMin(m_h),
m_hMax(m_h), m_kMin(m_k), m_kMax(m_k), m_lMin(m_l), m_lMax(m_l) {}
/// Return an iterator that can move from min to max, with current = min
HKLGenerator::const_iterator::const_iterator(const V3D &hklMin,
const V3D &hklMax)
: m_h(static_cast<int>(hklMin.X())), m_k(static_cast<int>(hklMin.Y())),
m_l(static_cast<int>(hklMin.Z())), m_hMin(m_h),
m_l(static_cast<int>(hklMin.Z())), m_hkl(hklMin), m_hMin(m_h),
m_hMax(static_cast<int>(hklMax.X())), m_kMin(m_k),
m_kMax(static_cast<int>(hklMax.Y())), m_lMin(m_l),
m_lMax(static_cast<int>(hklMax.Z())) {}
HKLGenerator::const_iterator::const_iterator(const V3D &hklMin,
const V3D &hklMax,
const V3D &current)
: m_h(static_cast<int>(current.X())), m_k(static_cast<int>(current.Y())),
m_l(static_cast<int>(current.Z())), m_hMin(static_cast<int>(hklMin.X())),
m_hMax(static_cast<int>(hklMax.X())),
m_kMin(static_cast<int>(hklMin.Y())),
m_kMax(static_cast<int>(hklMax.Y())),
m_lMin(static_cast<int>(hklMin.Z())),
m_lMax(static_cast<int>(hklMax.Z())) {}
/// Increments HKL, l moves fastest, h moves slowest, wrapping around at the max
/// for each index.
void HKLGenerator::const_iterator::increment() {
if (++m_l > m_lMax) {
++m_l;
if (m_l > m_lMax) {
m_l = m_lMin;
if (++m_k > m_kMax) {
++m_k;
if (m_k > m_kMax) {
m_k = m_kMin;
++m_h;
}
}
m_hkl.setX(m_h);
m_hkl.setY(m_k);
m_hkl.setZ(m_l);
m_hkl = V3D(m_h, m_k, m_l);
}
} // namespace Geometry
......
......@@ -38,64 +38,160 @@ public:
TS_ASSERT_EQUALS(gen.size(), 144);
}
void testSpeed() {
Timer t;
size_t N = 10;
void test_HKLGeneratorReturnsCorrectSizeOne() {
HKLGenerator gen(V3D(-2, -1, -5), V3D(-2, -1, -5));
TS_ASSERT_EQUALS(gen.size(), 1);
}
void test_beginIterator() {
HKLGenerator gen(V3D(-2, -2, -3), V3D(1, 1, 0));
HKLGenerator::const_iterator it = gen.begin();
TS_ASSERT_EQUALS(*it, V3D(-2, -2, -3));
}
void test_endIterator() {
HKLGenerator gen(V3D(-2, -2, -3), V3D(1, 1, 0));
HKLGenerator::const_iterator it = gen.end();
TS_ASSERT_EQUALS(*it, V3D(2, -2, -3));
}
void test_iterator_dereference() {
HKLGenerator::const_iterator it;
TS_ASSERT_EQUALS(*it, V3D(0, 0, 0));
}
SpaceGroup_const_sptr sg =
SpaceGroupFactory::Instance().createSpaceGroup("C m c m");
void test_comparison() {
HKLGenerator::const_iterator it1(V3D(-2, -3, 1));
HKLGenerator::const_iterator it2(V3D(-2, -3, 1));
HKLGenerator::const_iterator it3(V3D(-2, -3, 0));
HKLFilterSpaceGroup filter(sg);
HKLGenerator gen(30, 30, 30);
TS_ASSERT_EQUALS(it1, it2);
TS_ASSERT_DIFFERS(it1, it3);
}
void test_iterator_increment() {
HKLGenerator::const_iterator defaultIterator;
TS_ASSERT_THROWS_NOTHING(++defaultIterator);
}
void test_iterator_dereference_range() {
HKLGenerator::const_iterator it(V3D(-1, -1, -1), V3D(1, 1, 1));
TS_ASSERT_EQUALS(*it, V3D(-1, -1, -1));
++it;
TS_ASSERT_EQUALS(*it, V3D(-1, -1, 0));
}
void test_hkl_range() {
HKLGenerator::const_iterator it(V3D(-1, -1, -1), V3D(1, 1, 1));
HKLGenerator::const_iterator end(V3D(2, -1, -1));
std::vector<V3D> hkls;
std::copy(it, end, std::back_inserter(hkls));
TS_ASSERT_EQUALS(hkls.size(), 27);
TS_ASSERT_EQUALS(hkls[0], V3D(-1, -1, -1));
TS_ASSERT_EQUALS(hkls[1], V3D(-1, -1, 0));
TS_ASSERT_EQUALS(hkls[2], V3D(-1, -1, 1));
TS_ASSERT_EQUALS(hkls[3], V3D(-1, 0, -1));
TS_ASSERT_EQUALS(hkls[26], V3D(1, 1, 1));
TS_ASSERT_EQUALS(std::distance(it, end), 27);
}
void test_speed() {
size_t N = 1000;
Timer t;
for (size_t i = 0; i < N; ++i) {
HKLGenerator gen(20, 20, 20);
std::vector<V3D> hkls;
hkls.reserve(gen.size());
auto end = gen.end();
for (auto hkl = gen.begin(); hkl != end; ++hkl) {
if (filter(*hkl)) {
hkls.push_back(*hkl);
}
for (auto it = gen.begin(), end = gen.end(); it != end; ++it) {
hkls.push_back(*it);
}
// std::copy_if(gen.begin(), gen.end(), std::back_inserter(hkls),
// filter);
TS_ASSERT_LESS_THAN(0, hkls.size());
// TS_ASSERT_EQUALS(hkls.size(), 61 * 61 * 61);
TS_ASSERT_EQUALS(hkls.size(), 41 * 41 * 41);
}
std::cout << t.elapsed() / N << std::endl;
float e = t.elapsed();
std::cout << "T: " << e / static_cast<float>(N) << std::endl;
}
void testSpeed_old() {
void test_speed_filter_constructor() {
size_t N = 1000;
Timer t;
size_t N = 10;
for (size_t i = 0; i < N; ++i) {
HKLGenerator gen(20, 20, 20);
SpaceGroup_const_sptr sg =
SpaceGroupFactory::Instance().createSpaceGroup("C m c m");
std::vector<V3D> hkls;
hkls.reserve(gen.size());
std::copy(gen.begin(), gen.end(), std::back_inserter(hkls));
// TS_ASSERT_EQUALS(hkls.size(), gen.size());
TS_ASSERT_EQUALS(hkls.size(), 41 * 41 * 41);
}
float e = t.elapsed();
std::cout << "T: " << e / static_cast<float>(N) << std::endl;
}
void test_speed_filter() {
size_t N = 1000;
Timer t;
for (size_t i = 0; i < N; ++i) {
double dMin = 0.5;
UnitCell cell(8, 8, 8);
HKLGenerator gen(cell, dMin);
std::vector<V3D> hkls;
hkls.reserve(61 * 61 * 61);
for (int h = -30; h <= 30; ++h) {
for (int k = -30; k <= 30; ++k) {
for (int l = -30; l <= 30; ++l) {
V3D hkl(h, k, l);
if (sg->isAllowedReflection(hkl)) {
hkls.push_back(hkl);
}
}
hkls.reserve(gen.size());
HKLFilterDRange dFilter(cell, dMin, 200.0);
for (auto it = gen.begin(); it != gen.end(); ++it) {
if (dFilter.isAllowed(*it)) {
hkls.push_back(*it);
}
}
TS_ASSERT_LESS_THAN(0, hkls.size());
// TS_ASSERT_EQUALS(hkls.size(), 61 * 61 * 61);
// TS_ASSERT_EQUALS(hkls.size(), gen.size());
TS_ASSERT_DIFFERS(hkls.size(), 0);
}
float e = t.elapsed();
std::cout << "T: " << e / static_cast<float>(N) << std::endl;
}
void test_speed_filter_stl() {
size_t N = 1000;
Timer t;
for (size_t i = 0; i < N; ++i) {
double dMin = 0.5;
UnitCell cell(8, 8, 8);
HKLGenerator gen(cell, dMin);
std::vector<V3D> hkls;
hkls.reserve(gen.size());
HKLFilter_const_sptr dFilter = boost::make_shared<HKLFilterDRange>(cell, dMin, 200.0);
std::remove_copy_if(gen.begin(), gen.end(), std::back_inserter(hkls), HKLFilterProxy(~dFilter));
// TS_ASSERT_EQUALS(hkls.size(), gen.size());
TS_ASSERT_DIFFERS(hkls.size(), 0);
}
std::cout << t.elapsed() / N << std::endl;
float e = t.elapsed();
std::cout << "T: " << e / static_cast<float>(N) << std::endl;
}
};
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment