IntegrateQLabEvents.h

Go to the documentation of this file.

// Mantid Repository : https://github.com/mantidproject/mantid
//
// Copyright © 2012 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 +
#pragma once
 
#include "MantidDataObjects/EventWorkspace.h"
#include "MantidDataObjects/Peak.h"
#include "MantidDataObjects/PeakShapeEllipsoid.h"
#include "MantidKernel/Matrix.h"
#include "MantidKernel/V3D.h"
#include "MantidMDAlgorithms/DllConfig.h"
 
#include <memory>
#include <tuple>
#include <unordered_map>
#include <vector>
 
namespace Mantid {
 
namespace Geometry {
class PeakShape;
}
 
namespace MDAlgorithms {
using Mantid::DataObjects::PeakShapeEllipsoid_const_sptr;
using Mantid::Geometry::PeakShape_const_sptr;
using Mantid::Kernel::V3D;
 
struct CellCoords {
  int64_t a;
  int64_t b;
  int64_t c;
 
  CellCoords(const V3D &q, const double cellSize)
      : a(static_cast<int64_t>(q[0] / cellSize)), b(static_cast<int64_t>(q[1] / cellSize)),
        c(static_cast<int64_t>(q[2] / cellSize)) {}
 
  bool isOrigin() { return !(a || b || c); }
 
  int64_t getHash() { return 1000000000000 * a + 100000000 * b + 10000 * c; }
 
  std::vector<int64_t> nearbyCellHashes() {
    std::vector<int64_t> neighbors;
    for (int64_t ia = a - 1; ia <= a + 1; ia++)
      for (int64_t ib = b - 1; ib <= b + 1; ib++)
        for (int64_t ic = c - 1; ic <= c + 1; ic++) {
          int64_t key = 1000000000000 * ia + 100000000 * ib + 10000 * ic;
          neighbors.emplace_back(key);
        }
    return neighbors;
  }
};
 
// [(weight, error), Q-vector], trimmed-down info for an event
using SlimEvent = std::pair<std::pair<double, double>, V3D>;
using SlimEvents = std::vector<SlimEvent>;
 
// A cell in partitioned QLab space containing one peak
struct OccupiedCell {
  size_t peakIndex;  // index of the peak within this cell
  V3D peakQ;         // QLab vector of the peak within this cell
  SlimEvents events; // events potentially closer than m_radius to the peak
};
 
class MANTID_MDALGORITHMS_DLL IntegrateQLabEvents {
public:
  IntegrateQLabEvents(const SlimEvents &peak_q_list, double radius,
                      const bool useOnePercentBackgroundCorrection = true);
 
  static bool isOrigin(const V3D &q, const double &cellSize);
 
  void setRadius(const double &radius);
 
  void addEvents(SlimEvents const &event_qs);
 
  PeakShape_const_sptr ellipseIntegrateEvents(const std::vector<V3D> &E1Vec, V3D const &peak_q, bool specify_size,
                                              double peak_radius, double back_inner_radius, double back_outer_radius,
                                              std::vector<double> &axes_radii, double &inti, double &sigi,
                                              std::pair<double, double> &backi);
 
  void populateCellsWithPeaks();
 
private:
  static std::pair<double, double> numInEllipsoid(SlimEvents const &events, std::vector<V3D> const &directions,
                                                  std::vector<double> const &sizes);
 
  static std::pair<double, double> numInEllipsoidBkg(SlimEvents const &events, std::vector<V3D> const &directions,
                                                     std::vector<double> const &sizes,
                                                     std::vector<double> const &sizesIn,
                                                     const bool useOnePercentBackgroundCorrection);
 
  static void makeCovarianceMatrix(SlimEvents const &events, Kernel::DblMatrix &matrix, double radius);
 
  static void getEigenVectors(Kernel::DblMatrix const &cov_matrix, std::vector<V3D> &eigen_vectors,
                              std::vector<double> &eigen_values);
 
  void addEvent(const SlimEvent event);
 
  PeakShapeEllipsoid_const_sptr ellipseIntegrateEvents(const std::vector<V3D> &E1Vec, V3D const &peak_q,
                                                       SlimEvents const &ev_list, std::vector<V3D> const &directions,
                                                       std::vector<double> const &sigmas, bool specify_size,
                                                       double peak_radius, double back_inner_radius,
                                                       double back_outer_radius, std::vector<double> &axes_radii,
                                                       double &inti, double &sigi, std::pair<double, double> &backi);
 
  double detectorQ(const std::vector<V3D> &E1Vec, const V3D &QLabFrame, const std::vector<double> &r);
 
  // Private data members
  double m_radius; // size of sphere to use for events around a peak
  const bool m_useOnePercentBackgroundCorrection;
  double m_cellSize;
  std::unordered_map<size_t, OccupiedCell> m_cellsWithPeaks;
  std::unordered_map<size_t, SlimEvents> m_cellsWithEvents;
};
 
} // namespace MDAlgorithms
 
} // namespace Mantid