CountReflections.cpp

Go to the documentation of this file.

// Mantid Repository : https://github.com/mantidproject/mantid
//
// Copyright © 2018 ISIS Rutherford Appleton Laboratory UKRI,
//   NScD Oak Ridge National Laboratory, European Spallation Source,
//   Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
// SPDX - License - Identifier: GPL - 3.0 +
#include "MantidCrystal/CountReflections.h"
#include "MantidCrystal/PeakStatisticsTools.h"
 
#include "MantidAPI/Sample.h"
#include "MantidAPI/WorkspaceProperty.h"
 
#include "MantidDataObjects/PeaksWorkspace.h"
 
#include "MantidGeometry/Crystal/OrientedLattice.h"
#include "MantidGeometry/Crystal/PointGroupFactory.h"
#include "MantidGeometry/Crystal/ReflectionCondition.h"
 
#include "MantidKernel/ListValidator.h"
 
namespace Mantid::Crystal {
 
using Mantid::Kernel::Direction;
 
using namespace Mantid::API;
using namespace Mantid::DataObjects;
using namespace Mantid::Geometry;
using namespace Mantid::Kernel;
 
// Register the algorithm into the AlgorithmFactory
DECLARE_ALGORITHM(CountReflections)
 
//----------------------------------------------------------------------------------------------
 
 
const std::string CountReflections::name() const { return "CountReflections"; }
 
int CountReflections::version() const { return 1; }
 
const std::string CountReflections::category() const { return "Crystal\\Peaks"; }
 
const std::string CountReflections::summary() const {
  return "Calculates statistics for a PeaksWorkspace based on symmetry and "
         "counting reflections.";
}
 
//----------------------------------------------------------------------------------------------
void CountReflections::init() {
  declareProperty(std::make_unique<WorkspaceProperty<PeaksWorkspace>>("InputWorkspace", "", Direction::Input),
                  "A workspace with peaks to calculate statistics for. Sample "
                  "with valid UB-matrix is required.");
 
  auto centeringSymbols = getAllReflectionConditionSymbols();
  declareProperty("LatticeCentering", centeringSymbols[0], std::make_shared<StringListValidator>(centeringSymbols),
                  "Lattice centering of the cell.");
 
  auto pointGroups = PointGroupFactory::Instance().getAllPointGroupSymbols();
  declareProperty("PointGroup", "1", std::make_shared<StringListValidator>(pointGroups),
                  "Point group symmetry for completeness and redundancy calculations.");
 
  declareProperty(std::make_unique<PropertyWithValue<double>>("MinDSpacing", 1.0, Direction::Input),
                  "Minimum d-spacing for completeness calculation.");
 
  declareProperty(std::make_unique<PropertyWithValue<double>>("MaxDSpacing", 100.0, Direction::Input),
                  "Maximum d-spacing for completeness calculation.");
 
  declareProperty(std::make_unique<PropertyWithValue<int>>("UniqueReflections", 0, Direction::Output),
                  "Number of unique reflections in data set.");
 
  declareProperty(std::make_unique<PropertyWithValue<double>>("Completeness", 0.0, Direction::Output),
                  "Completeness of the data set as a fraction between 0 and 1.");
 
  declareProperty(std::make_unique<PropertyWithValue<double>>("Redundancy", 0.0, Direction::Output),
                  "Average redundancy in data set, depending on point group.");
 
  declareProperty(std::make_unique<PropertyWithValue<double>>("MultiplyObserved", 0.0, Direction::Output),
                  "Fraction of reflections with more than one observation.");
 
  declareProperty(std::make_unique<WorkspaceProperty<PeaksWorkspace>>("MissingReflectionsWorkspace", "",
                                                                      Direction::Output, PropertyMode::Optional),
                  "Reflections in specified d-range that are missing in input workspace.");
}
 
//----------------------------------------------------------------------------------------------
void CountReflections::exec() {
  double dMin = getProperty("MinDSpacing");
  double dMax = getProperty("MaxDSpacing");
 
  PointGroup_sptr pointGroup = PointGroupFactory::Instance().createPointGroup(getProperty("PointGroup"));
 
  ReflectionCondition_sptr centering = getReflectionConditionBySymbol(getProperty("LatticeCentering"));
 
  PeaksWorkspace_sptr inputPeaksWorkspace = getProperty("InputWorkspace");
 
  UnitCell cell = inputPeaksWorkspace->sample().getOrientedLattice();
 
  PeakStatisticsTools::UniqueReflectionCollection reflections(cell, std::make_pair(dMin, dMax), pointGroup, centering);
 
  auto peaks = inputPeaksWorkspace->getPeaks();
  reflections.addObservations(peaks);
 
  auto possibleUniqueReflections = static_cast<double>(reflections.getUniqueReflectionCount());
 
  size_t observedUniqueReflections = reflections.getObservedUniqueReflectionCount();
 
  auto observedUniqueReflectionsD = static_cast<double>(observedUniqueReflections);
 
  size_t totalReflections = reflections.getObservedReflectionCount();
 
  if (peaks.size() > totalReflections) {
    g_log.information() << "There are " << (peaks.size() - totalReflections)
                        << " peaks in the input workspace that fall outside "
                           "the resolution limit and are not considered for "
                           "the calculations."
                        << std::endl;
  }
 
  auto multiplyObservedReflections = static_cast<double>(reflections.getObservedUniqueReflectionCount(1));
 
  setProperty("UniqueReflections", static_cast<int>(observedUniqueReflections));
  setProperty("Completeness", observedUniqueReflectionsD / possibleUniqueReflections);
  setProperty("Redundancy", static_cast<double>(totalReflections) / observedUniqueReflectionsD);
  setProperty("MultiplyObserved", multiplyObservedReflections / observedUniqueReflectionsD);
 
  PeaksWorkspace_sptr outputWorkspace = getPeaksWorkspace(inputPeaksWorkspace, reflections, pointGroup);
 
  if (outputWorkspace) {
    setProperty("MissingReflectionsWorkspace", outputWorkspace);
  }
}
 
PeaksWorkspace_sptr
CountReflections::getPeaksWorkspace(const PeaksWorkspace_sptr &templateWorkspace,
                                    const PeakStatisticsTools::UniqueReflectionCollection &reflections,
                                    const PointGroup_sptr &pointGroup) const {
  std::string outputWorkspaceName = getPropertyValue("MissingReflectionsWorkspace");
 
  if (outputWorkspaceName.empty()) {
    return PeaksWorkspace_sptr();
  }
 
  PeaksWorkspace_sptr rawOutputPeaksWorkspace = getProperty("MissingReflectionsWorkspace");
 
  PeaksWorkspace_sptr outputPeaksWorkspace = std::dynamic_pointer_cast<PeaksWorkspace>(rawOutputPeaksWorkspace);
 
  if (outputPeaksWorkspace != templateWorkspace) {
    outputPeaksWorkspace = templateWorkspace->clone();
  }
 
  const auto &missingPeaks = reflections.getUnobservedUniqueReflections();
 
  std::vector<Peak> peaks;
  peaks.reserve(missingPeaks.size() * pointGroup->order());
 
  for (const auto &reflection : missingPeaks) {
    auto hkls = pointGroup->getEquivalents(reflection);
 
    for (const auto &hkl : hkls) {
      Peak peak;
      peak.setHKL(hkl);
 
      peaks.emplace_back(peak);
    }
  }
 
  outputPeaksWorkspace->getPeaks().swap(peaks);
 
  return outputPeaksWorkspace;
}
 
} // namespace Mantid::Crystal