Mantid
Loading...
Searching...
No Matches
DiscusMultipleScatteringCorrection.h
Go to the documentation of this file.
1// Mantid Repository : https://github.com/mantidproject/mantid
2//
3// Copyright © 2020 ISIS Rutherford Appleton Laboratory UKRI,
4// NScD Oak Ridge National Laboratory, European Spallation Source,
5// Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
6// SPDX - License - Identifier: GPL - 3.0 +
7#pragma once
8
9//------------------------------------------------------------------------------
10// Includes
11//------------------------------------------------------------------------------
12#include "MantidAPI/Algorithm.h"
13#include "MantidAlgorithms/DllConfig.h"
20
21#include <boost/container/small_vector.hpp>
22
23namespace Mantid {
24namespace API {
25class Sample;
26}
27namespace Geometry {
28class Instrument;
29}
30
31namespace Algorithms {
32
33// define some simple classes to store 2D datasets instead of using MatrixWorkspace internally
34// This couples the algorithm more loosely to Mantid and avoids some complexity in choosing whether
35// to call readX, dataX etc
37 // separate vectors of X and Y rather than vector of pairs to mirror Histogram class and support edges\points
38 std::vector<double> X;
39 std::vector<double> Y;
41 DiscusData1D(std::vector<double> X, std::vector<double> Y) : X(std::move(X)), Y(std::move(Y)) {}
42};
43
45public:
46 DiscusData2D() : m_data(std::vector<DiscusData1D>{}), m_specAxis(nullptr){};
47 DiscusData2D(const std::vector<DiscusData1D> &data, const std::shared_ptr<std::vector<double>> &specAxis)
48 : m_data(data), m_specAxis(specAxis){};
49 std::unique_ptr<DiscusData2D> createCopy(bool clearY = false);
50 size_t getNumberHistograms() { return m_data.size(); }
51 DiscusData1D &histogram(const size_t i) { return m_data[i]; }
52 std::vector<DiscusData1D> &histograms() { return m_data; }
53 const std::vector<double> &getSpecAxisValues();
54
55private:
56 std::vector<DiscusData1D> m_data;
57 // optional spectrum axis
58 std::shared_ptr<std::vector<double>> m_specAxis;
59};
60
63 std::string_view materialName;
64 std::shared_ptr<DiscusData2D> SQ;
65 std::shared_ptr<DiscusData2D> logSQ{};
66 std::shared_ptr<DiscusData1D> QSQScaleFactor{};
67 std::shared_ptr<DiscusData2D> QSQ{};
68 std::shared_ptr<DiscusData2D> InvPOfQ{};
69 std::shared_ptr<int> scatterCount = std::make_shared<int>(0);
70};
71
78class MANTID_ALGORITHMS_DLL DiscusMultipleScatteringCorrection : public API::Algorithm {
79public:
80 // use small_vector to avoid performance hit from heap allocation of std::vector. Use size 5 in line with Track.h
81 using ComponentWorkspaceMappings = boost::container::small_vector<ComponentWorkspaceMapping, 5>;
84 const std::string name() const override { return "DiscusMultipleScatteringCorrection"; }
86 int version() const override { return 1; }
87 const std::vector<std::string> seeAlso() const override {
88 return {"MayersSampleCorrection", "CarpenterSampleCorrection", "VesuvioCalculateMS"};
89 }
91 const std::string category() const override { return "CorrectionFunctions"; }
93 const std::string summary() const override {
94 return "Calculates a multiple scattering correction using a Monte Carlo method";
95 }
96 const std::string alias() const override { return "Muscat"; }
97 bool checkGroups() override { return false; }
98
99protected:
100 virtual std::shared_ptr<SparseWorkspace> createSparseWorkspace(const API::MatrixWorkspace &modelWS,
101 const size_t nXPoints, const size_t rows,
102 const size_t columns);
103 virtual std::unique_ptr<InterpolationOption> createInterpolateOption();
104 double interpolateFlat(const DiscusData1D &histToInterpolate, double x);
105 std::tuple<double, int> sampleQW(const std::shared_ptr<DiscusData2D> &CumulativeProb, double x);
106 double interpolateSquareRoot(const DiscusData1D &histToInterpolate, double x);
107 double interpolateGaussian(const DiscusData1D &histToInterpolate, double x);
108 double Interpolate2D(const ComponentWorkspaceMapping &SQWSMapping, double q, double w);
109 void updateTrackDirection(Geometry::Track &track, const double cosT, const double phi);
110 void integrateCumulative(const DiscusData1D &h, const double xmin, const double xmax, std::vector<double> &resultX,
111 std::vector<double> &resultY, const bool returnCumulative);
113 void getXMinMax(const Mantid::API::MatrixWorkspace &ws, double &xmin, double &xmax) const;
114 void prepareSampleBeamGeometry(const API::MatrixWorkspace_sptr &inputWS);
115
116private:
117 void init() override;
118 void exec() override;
119 std::map<std::string, std::string> validateInputs() override;
120 API::MatrixWorkspace_sptr createOutputWorkspace(const API::MatrixWorkspace &inputWS) const;
121 std::tuple<double, double> new_vector(const Kernel::Material &material, double k, bool specialSingleScatterCalc);
122 std::vector<double> simulatePaths(const int nEvents, const int nScatters, Kernel::PseudoRandomNumberGenerator &rng,
123 ComponentWorkspaceMappings &componentWorkspaces, const double kinc,
124 const std::vector<double> &wValues, const Kernel::V3D &detPos,
125 bool specialSingleScatterCalc);
126 std::tuple<bool, std::vector<double>> scatter(const int nScatters, Kernel::PseudoRandomNumberGenerator &rng,
127 const ComponentWorkspaceMappings &componentWorkspaces,
128 const double kinc, const std::vector<double> &wValues,
129 const Kernel::V3D &detPos, bool specialSingleScatterCalc);
131 Geometry::Track generateInitialTrack(Kernel::PseudoRandomNumberGenerator &rng);
132 void inc_xyz(Geometry::Track &track, double vl);
133 const Geometry::IObject *updateWeightAndPosition(Geometry::Track &track, double &weight, const double k,
135 bool specialSingleScatterCalc,
136 const ComponentWorkspaceMappings &componentWorkspaces);
137 bool q_dir(Geometry::Track &track, const Geometry::IObject *shapePtr, const ComponentWorkspaceMappings &invPOfQs,
138 double &k, const double scatteringXSection, Kernel::PseudoRandomNumberGenerator &rng, double &weight);
139 void interpolateFromSparse(API::MatrixWorkspace &targetWS, const SparseWorkspace &sparseWS,
141 void correctForWorkspaceNameClash(std::string &wsName);
142 void setWorkspaceName(const API::MatrixWorkspace_sptr &ws, std::string wsName);
143 void createInvPOfQWorkspaces(ComponentWorkspaceMappings &matWSs, size_t nhists);
144 void convertToLogWorkspace(const std::shared_ptr<DiscusData2D> &SOfQ);
145 void calculateQSQIntegralAsFunctionOfK(ComponentWorkspaceMappings &matWSs, const std::vector<double> &specialKs);
146 void prepareCumulativeProbForQ(double kinc, const ComponentWorkspaceMappings &PInvOfQs);
147 void prepareQSQ(double kinc);
148 double getKf(const double deltaE, const double kinc);
149 std::tuple<double, double, int, double> sampleQWUniform(const std::vector<double> &wValues,
150 Kernel::PseudoRandomNumberGenerator &rng, const double kinc);
151 void prepareStructureFactors();
152 void convertWsBothAxesToPoints(API::MatrixWorkspace_sptr &ws);
153 std::tuple<double, double> getKinematicRange(double kf, double ki);
154 std::vector<std::tuple<double, int, double>> generateInputKOutputWList(const double efixed,
155 const std::vector<double> &xPoints);
156 std::tuple<std::vector<double>, std::vector<double>, std::vector<double>>
157 integrateQSQ(const std::shared_ptr<DiscusData2D> &QSQ, double kinc, const bool returnCumulative);
158 double getQSQIntegral(const DiscusData1D &QSQScaleFactor, double k);
159 const ComponentWorkspaceMapping *findMatchingComponent(const ComponentWorkspaceMappings &componentWorkspaces,
160 const Geometry::IObject *shapeObjectWithScatter);
161 void addWorkspaceToDiscus2DData(const Geometry::IObject_const_sptr &shape, const std::string_view &matName,
163 long long m_callsToInterceptSurface{0};
164 long long m_IkCalculations{0};
166 int m_maxScatterPtAttempts{};
167 std::shared_ptr<const DiscusData1D> m_sigmaSS; // scattering cross section as a function of k
168 // vectors of S(Q,w) and derived quantities. One entry for sample and each environment component
171 bool m_importanceSampling{};
172 Kernel::DeltaEMode::Type m_EMode{Kernel::DeltaEMode::Undefined};
173 bool m_simulateEnergiesIndependently{};
175 std::shared_ptr<const Geometry::ReferenceFrame> m_refframe;
176 const Geometry::SampleEnvironment *m_env{nullptr};
177 bool m_NormalizeSQ{};
179 std::unique_ptr<IBeamProfile> m_beamProfile;
180};
181} // namespace Algorithms
182} // namespace Mantid
Base class from which all concrete algorithm classes should be derived.
Definition: Algorithm.h:85
Base MatrixWorkspace Abstract Class.
DiscusData2D(const std::vector< DiscusData1D > &data, const std::shared_ptr< std::vector< double > > &specAxis)
std::shared_ptr< std::vector< double > > m_specAxis
std::unique_ptr< DiscusData2D > createCopy(bool clearY=false)
Calculates a multiple scattering correction Based on Muscat Fortran code provided by Spencer Howells.
const std::vector< std::string > seeAlso() const override
Function to return all of the seeAlso algorithms related to this algorithm.
const std::string alias() const override
function to return any aliases of the algorithm.
const std::string category() const override
Algorithm's category for identification.
const std::string summary() const override
Summary of algorithms purpose.
const std::string name() const override
Algorithm's name.
bool checkGroups() override
Check the input workspace properties for groups.
boost::container::small_vector< ComponentWorkspaceMapping, 5 > ComponentWorkspaceMappings
std::shared_ptr< const Geometry::ReferenceFrame > m_refframe
Class to provide a consistent interface to an interpolation option on algorithms.
Defines functions and utilities to create and deal with sparse instruments.
A simple structure that defines an axis-aligned cuboid shaped bounding box for a geometrical object.
Definition: BoundingBox.h:34
IObject : Interface for geometry objects.
Definition: IObject.h:41
Defines a single instance of a SampleEnvironment.
Defines a track as a start point and a direction.
Definition: Track.h:165
A material is defined as being composed of a given element, defined as a PhysicalConstants::NeutronAt...
Definition: Material.h:50
Defines a 1D pseudo-random number generator, i.e.
Class for 3D vectors.
Definition: V3D.h:34
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
std::shared_ptr< const IObject > IObject_const_sptr
Typdef for a shared pointer to a const object.
Definition: IObject.h:94
Helper class which provides the Collimation Length for SANS instruments.
STL namespace.
DiscusData1D(std::vector< double > X, std::vector< double > Y)
Type
Define the available energy transfer modes It is important to assign enums proper numbers,...
Definition: DeltaEMode.h:29