Mantid
Loading...
Searching...
No Matches
DoublePulseFit.cpp
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 +
9
16#include "MantidAPI/TableRow.h"
20
24
25#include <math.h>
26#include <memory>
27
29
31
33 const Mantid::API::MatrixWorkspace_sptr &workspace, const std::string &suffix) {
34 fit.setProperty("InputWorkspace" + suffix, workspace);
35
36 int workspaceIndex = fittingAlgorithm.getProperty("WorkspaceIndex" + suffix);
37 fit.setProperty("WorkspaceIndex" + suffix, workspaceIndex);
38
39 double startX = fittingAlgorithm.getProperty("StartX" + suffix);
40 double endX = fittingAlgorithm.getProperty("EndX" + suffix);
41 fit.setProperty("StartX" + suffix, startX);
42 fit.setProperty("EndX" + suffix, endX);
43
44 std::vector<double> exclude = fittingAlgorithm.getProperty("Exclude" + suffix);
45 fit.setProperty("Exclude" + suffix, exclude);
46}
47
49 const std::vector<Mantid::API::MatrixWorkspace_sptr> &workspaces) {
50 setMultiDataProperties(fittingAlgorithm, fit, workspaces[0], "");
51
52 for (auto i = 1u; i < workspaces.size(); ++i)
53 setMultiDataProperties(fittingAlgorithm, fit, workspaces[i], "_" + std::to_string(i));
54}
55
56Mantid::API::IFunction_sptr getDoublePulseFunction(std::shared_ptr<const API::IFunction> const &function, double offset,
57 double firstPulseWeight, double secondPulseWeight) {
58 auto clonedFunction = function->clone();
59
60 auto convolution = std::make_shared<Mantid::CurveFitting::Functions::Convolution>();
61 auto delta1 = std::make_shared<Mantid::CurveFitting::Functions::DeltaFunction>();
62 auto delta2 = std::make_shared<Mantid::CurveFitting::Functions::DeltaFunction>();
63 auto deltaComposite = std::make_shared<Mantid::API::CompositeFunction>();
64
65 convolution->setAttributeValue("FixResolution", false);
66 delta1->setParameter("Centre", -0.5 * offset);
67 delta1->setParameter("Height", firstPulseWeight);
68 delta1->fixAll();
69 delta2->setParameter("Centre", 0.5 * offset);
70 delta2->setParameter("Height", secondPulseWeight);
71 delta2->fixAll();
72 deltaComposite->addFunction(delta1);
73 deltaComposite->addFunction(delta2);
74
75 convolution->addFunction(clonedFunction);
76 convolution->addFunction(deltaComposite);
77
78 return convolution;
79}
80
82getDoublePulseMultiDomainFunction(std::shared_ptr<const API::MultiDomainFunction> const &initialFunction, double offset,
83 double firstPulseWeight, double secondPulseWeight) {
84 auto doublePulseFunc = std::make_shared<API::MultiDomainFunction>();
85 for (size_t domain = 0; domain < initialFunction->getNumberDomains(); domain++) {
86 auto twoPulseInnerFunction =
87 getDoublePulseFunction(initialFunction->getFunction(domain), offset, firstPulseWeight, secondPulseWeight);
88 doublePulseFunc->addFunction(twoPulseInnerFunction);
89 doublePulseFunc->setDomainIndex(domain, domain);
90 }
91 return doublePulseFunc;
92}
93
95extractInnerFunction(std::shared_ptr<const Mantid::CurveFitting::Functions::Convolution> const &convolutionFunction) {
96 return convolutionFunction->getFunction(0);
97}
98
100extractInnerFunction(std::shared_ptr<const API::MultiDomainFunction> const &doublePulseFunc) {
101 auto extractedFunction = std::make_shared<API::MultiDomainFunction>();
102 for (size_t domain = 0; domain < doublePulseFunc->getNumberDomains(); domain++) {
103 auto convFunction = std::dynamic_pointer_cast<Convolution>(doublePulseFunc->getFunction(domain));
104 if (convFunction) {
105 extractedFunction->addFunction(extractInnerFunction(convFunction));
106 extractedFunction->setDomainIndex(domain, domain);
107 } else {
108 throw std::runtime_error("Cannot extract from non convolution function. DoublePulseFit.cpp");
109 }
110 }
111 return extractedFunction;
112}
113
114// Register the class into the algorithm factory
116
117
119
124 getPointerToProperty("Ties")->setDocumentation("Math expressions defining ties between parameters of "
125 "the fitting function.");
126 declareProperty("Constraints", "", Kernel::Direction::Input);
127 getPointerToProperty("Constraints")->setDocumentation("List of constraints");
128 auto mustBePositive = std::make_shared<Kernel::BoundedValidator<int>>();
129 mustBePositive->setLower(0);
130 declareProperty("MaxIterations", 500, mustBePositive->clone(),
131 "Stop after this number of iterations if a good fit is not found");
132 declareProperty("OutputStatus", "", Kernel::Direction::Output);
133 getPointerToProperty("OutputStatus")->setDocumentation("Whether the fit was successful");
134 declareProperty("OutputChi2overDoF", 0.0, "Returns the goodness of the fit", Kernel::Direction::Output);
135
136 std::vector<std::string> minimizerOptions = API::FuncMinimizerFactory::Instance().getKeys();
137 Kernel::IValidator_sptr minimizerValidator = std::make_shared<Kernel::StartsWithValidator>(minimizerOptions);
138
139 declareProperty("Minimizer", "Levenberg-Marquardt", minimizerValidator, "Minimizer to use for fitting.");
140
141 std::vector<std::string> costFuncOptions = API::CostFunctionFactory::Instance().getKeys();
142 // select only CostFuncFitting variety
143 for (auto &costFuncOption : costFuncOptions) {
144 auto costFunc = std::dynamic_pointer_cast<CostFunctions::CostFuncFitting>(
146 if (!costFunc) {
147 costFuncOption = "";
148 }
149 }
150 Kernel::IValidator_sptr costFuncValidator = std::make_shared<Kernel::ListValidator<std::string>>(costFuncOptions);
151 declareProperty("CostFunction", "Least squares", costFuncValidator,
152 "The cost function to be used for the fit, default is Least squares", Kernel::Direction::InOut);
153 declareProperty("CreateOutput", false,
154 "Set to true to create output workspaces with the results of the fit"
155 "(default is false).");
156 declareProperty("Output", "",
157 "A base name for the output workspaces (if not "
158 "given default names will be created). The "
159 "default is to use the name of the original data workspace as prefix "
160 "followed by suffixes _Workspace, _Parameters, etc.");
161 declareProperty("CalcErrors", false,
162 "Set to true to calcuate errors when output isn't created "
163 "(default is false).");
164 declareProperty("OutputCompositeMembers", false,
165 "If true and CreateOutput is true then the value of each "
166 "member of a Composite Function is also output.");
167 declareProperty(std::make_unique<Kernel::PropertyWithValue<bool>>("ConvolveMembers", false),
168 "If true and OutputCompositeMembers is true members of any "
169 "Convolution are output convolved\n"
170 "with corresponding resolution");
171 declareProperty("OutputParametersOnly", false,
172 "Set to true to output only the parameters and not "
173 "workspace(s) with the calculated values\n"
174 "(default is false, ignored if CreateOutput is false and "
175 "Output is an empty string).");
176 declareProperty("PulseOffset", 0.0, "The time offset between the two pulses");
177 declareProperty("FirstPulseWeight", 0.5, "Weighting of first pulse.");
178 declareProperty("SecondPulseWeight", 0.5, "Weighting of second pulse.");
179}
180
184
186 Mantid::API::IFunction_sptr function = getProperty("Function");
187 Mantid::API::IFunction_sptr doublePulseFunction;
188 auto pulseOffset = getProperty("PulseOffset");
189 auto firstPulseWeight = getProperty("FirstPulseWeight");
190 auto secondPulseWeight = getProperty("SecondPulseWeight");
191 if (auto multiDomainFunction = std::dynamic_pointer_cast<Mantid::API::MultiDomainFunction>(function)) {
192 doublePulseFunction =
193 getDoublePulseMultiDomainFunction(multiDomainFunction, pulseOffset, firstPulseWeight, secondPulseWeight);
194 } else {
195 doublePulseFunction = getDoublePulseFunction(function, pulseOffset, firstPulseWeight, secondPulseWeight);
196 }
197
198 runFitAlgorith(fitAlg, doublePulseFunction, getProperty("MaxIterations"));
199
200 createOutput(fitAlg, doublePulseFunction);
201}
202
210 std::string baseName = getProperty("Output");
211 if (baseName.empty()) {
212 API::Workspace_const_sptr ws = getProperty("InputWorkspace");
213 baseName = ws->getName();
214 }
215
216 if (m_makeOutput) {
217 declareProperty(std::make_unique<API::WorkspaceProperty<API::ITableWorkspace>>("OutputNormalisedCovarianceMatrix",
219 "The name of the TableWorkspace in which to store the final "
220 "covariance matrix");
221 setPropertyValue("OutputNormalisedCovarianceMatrix", baseName + "_NormalisedCovarianceMatrix");
222
223 declareProperty(std::make_unique<API::WorkspaceProperty<API::ITableWorkspace>>("OutputParameters", "",
225 "The name of the TableWorkspace in which to store the "
226 "final fit parameters");
227
228 setPropertyValue("OutputParameters", baseName + "_Parameters");
229
230 if (m_outputFitData) {
231 if (m_multiDomain) {
233 "OutputWorkspace", "", Kernel::Direction::Output),
234 "Name of the output Workspace holding resulting fitted "
235 "spectrum");
236 setPropertyValue("OutputWorkspace", baseName + "_Workspace");
237 } else {
239 "OutputWorkspace", "", Kernel::Direction::Output),
240 "Name of the output Workspace holding resulting simulated "
241 "spectrum");
242 setPropertyValue("OutputWorkspace", baseName + "_Workspace");
243 }
244 }
245 }
246}
247
248std::vector<Mantid::API::MatrixWorkspace_sptr> DoublePulseFit::getWorkspaces() const {
249 std::vector<Mantid::API::MatrixWorkspace_sptr> workspaces;
250 workspaces.reserve(m_workspacePropertyNames.size());
251 for (const auto &propertyName : m_workspacePropertyNames) {
253 workspaces.emplace_back(std::dynamic_pointer_cast<Mantid::API::MatrixWorkspace>(workspace));
254 }
255 return workspaces;
256}
257
259 std::string baseName = getProperty("Output");
260 bool makeOutput = getProperty("CreateOutput");
261 m_makeOutput = makeOutput || !baseName.empty();
262 m_outputFitData = !getProperty("OutputParametersOnly");
263 Mantid::API::IFunction_sptr function = getProperty("Function");
264 m_multiDomain = function->getNumberDomains() > 1;
265}
266
268 const Mantid::API::IFunction_sptr &function, int maxIterations) {
269 const bool convolveMembers = getProperty("ConvolveMembers");
270 const bool ignoreInvalidData = getProperty("IgnoreInvalidData");
271 const bool calcErrors = getProperty("CalcErrors");
272 const bool createOutput = getProperty("CreateOutput");
273 const bool outputCompositeMembers = getProperty("OutputCompositeMembers");
274 const bool outputParametersOnly = getProperty("OutputParametersOnly");
275
276 fitAlg->setLogging(true);
277 fitAlg->setProperty("Function", function);
278 setMultiDataProperties(*this, *fitAlg, getWorkspaces());
279 fitAlg->setProperty("IgnoreInvalidData", ignoreInvalidData);
280 fitAlg->setProperty("DomainType", getPropertyValue("DomainType"));
281 fitAlg->setProperty("EvaluationType", getPropertyValue("EvaluationType"));
282 fitAlg->setPropertyValue("PeakRadius", getPropertyValue("PeakRadius"));
283 fitAlg->setProperty("Ties", getPropertyValue("Ties"));
284 fitAlg->setProperty("Constraints", getPropertyValue("Constraints"));
285 fitAlg->setProperty("MaxIterations", maxIterations);
286 fitAlg->setProperty("Minimizer", getPropertyValue("Minimizer"));
287 fitAlg->setProperty("CostFunction", getPropertyValue("CostFunction"));
288 fitAlg->setProperty("CalcErrors", calcErrors);
289 fitAlg->setProperty("OutputCompositeMembers", outputCompositeMembers);
290 fitAlg->setProperty("ConvolveMembers", convolveMembers);
291 fitAlg->setProperty("CreateOutput", createOutput);
292 fitAlg->setProperty("OutputParametersOnly", outputParametersOnly);
293 fitAlg->setProperty("Output", getPropertyValue("Output"));
294 fitAlg->executeAsChildAlg();
295}
296
298 const Mantid::API::IFunction_sptr &function) {
299 Mantid::API::IFunction_sptr extractedFunction;
300 if (auto convFunction = std::dynamic_pointer_cast<Mantid::CurveFitting::Functions::Convolution>(function)) {
301 extractedFunction = extractInnerFunction(convFunction);
302 } else if (auto multiDomainFunction = std::dynamic_pointer_cast<Mantid::API::MultiDomainFunction>(function)) {
303 extractedFunction = extractInnerFunction(multiDomainFunction);
304 } else {
305 throw std::runtime_error("Incompatible function form. DoublePulseFit.cpp");
306 }
307 double outputChi2overDoF = fitAlg->getProperty("OutputChi2overDoF");
308 std::string outputStatus = fitAlg->getProperty("OutputStatus");
309 setProperty("OutputStatus", outputStatus);
310 setProperty("OutputChi2overDoF", outputChi2overDoF);
311 if (m_makeOutput) {
312 Mantid::API::ITableWorkspace_sptr covarianceMatrix = fitAlg->getProperty("OutputNormalisedCovarianceMatrix");
313 setProperty("OutputNormalisedCovarianceMatrix", covarianceMatrix);
314 Mantid::API::ITableWorkspace_sptr parametersTable = fitAlg->getProperty("OutputParameters");
315 setProperty("OutputParameters", parametersTable);
316 if (m_outputFitData) {
317 if (m_multiDomain) {
318 Mantid::API::WorkspaceGroup_sptr outputWorkspace = fitAlg->getProperty("OutputWorkspace");
319 setProperty("OutputWorkspace", outputWorkspace);
320 } else {
321 Mantid::API::MatrixWorkspace_sptr outputWorkspace = fitAlg->getProperty("OutputWorkspace");
322 setProperty("OutputWorkspace", outputWorkspace);
323 }
324 }
325 }
326 // // Fit has the function property as an InOut object and modifies it in
327 // place.
328 // // To replicate this behaviour we reset the relevant pointer to the new
329 // // function.
330 // Mantid::API::IFunction_sptr inputFunction = getProperty("Function");
331 // inputFunction.reset(extractedFunction.get());
332 // extractedFunction.reset();
333 setProperty("Function", extractedFunction);
334}
335
336} // namespace Mantid::CurveFitting::Algorithms
#define DECLARE_ALGORITHM(classname)
Definition: Algorithm.h:576
IPeaksWorkspace_sptr workspace
Definition: IndexPeaks.cpp:114
void declareProperty(std::unique_ptr< Kernel::Property > p, const std::string &doc="") override
Add a property to the list of managed properties.
Definition: Algorithm.cpp:1913
Kernel::Property * getPointerToProperty(const std::string &name) const override
Get a property by name.
Definition: Algorithm.cpp:2033
std::string getPropertyValue(const std::string &name) const override
Get the value of a property as a string.
Definition: Algorithm.cpp:2026
TypedValue getProperty(const std::string &name) const override
Get the value of a property.
Definition: Algorithm.cpp:2076
virtual std::shared_ptr< Algorithm > createChildAlgorithm(const std::string &name, const double startProgress=-1., const double endProgress=-1., const bool enableLogging=true, const int &version=-1)
Create a Child Algorithm.
Definition: Algorithm.cpp:842
void setPropertyValue(const std::string &name, const std::string &value) override
Set the value of a property by string N.B.
Definition: Algorithm.cpp:1975
IAlgorithm is the interface implemented by the Algorithm base class.
Definition: IAlgorithm.h:45
A property class for workspaces.
A function to fit muon data from a double pulse source.
void createOutput(const Mantid::API::IAlgorithm_sptr &fitAlg, const Mantid::API::IFunction_sptr &function)
void runFitAlgorith(const Mantid::API::IAlgorithm_sptr &fitAlgorithm, const Mantid::API::IFunction_sptr &function, int maxIterations)
std::vector< Mantid::API::MatrixWorkspace_sptr > getWorkspaces() const
void initConcrete() override
Initialisation method.
void execConcrete() override
Child classes implement the algorithm logic here.
void declareAdditionalProperties()
This function declares the additional output properties which are not present in init.
A base class for fitting algorithms.
std::vector< std::string > m_workspacePropertyNames
IPropertyManager * setProperty(const std::string &name, const T &value)
Templated method to set the value of a PropertyWithValue.
virtual TypedValue getProperty(const std::string &name) const =0
Get the value of a property.
The concrete, templated class for properties.
void setDocumentation(const std::string &documentation)
Sets the user level description of the property.
Definition: Property.cpp:134
static T & Instance()
Return a reference to the Singleton instance, creating it if it does not already exist Creation is do...
std::shared_ptr< IAlgorithm > IAlgorithm_sptr
shared pointer to Mantid::API::IAlgorithm
std::shared_ptr< WorkspaceGroup > WorkspaceGroup_sptr
shared pointer to Mantid::API::WorkspaceGroup
std::shared_ptr< ITableWorkspace > ITableWorkspace_sptr
shared pointer to Mantid::API::ITableWorkspace
std::shared_ptr< Workspace > Workspace_sptr
shared pointer to Mantid::API::Workspace
Definition: Workspace_fwd.h:20
std::shared_ptr< const Workspace > Workspace_const_sptr
shared pointer to Mantid::API::Workspace (const version)
Definition: Workspace_fwd.h:22
std::shared_ptr< IFunction > IFunction_sptr
shared pointer to the function base class
Definition: IFunction.h:732
std::shared_ptr< Algorithm > Algorithm_sptr
Typedef for a shared pointer to an Algorithm.
Definition: Algorithm.h:61
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
void setMultiDataProperties(const Mantid::API::IAlgorithm &fittingAlgorithm, Mantid::API::IAlgorithm &fit, const Mantid::API::MatrixWorkspace_sptr &workspace, const std::string &suffix)
MANTID_CURVEFITTING_DLL Mantid::API::IFunction_sptr getDoublePulseFunction(std::shared_ptr< const API::IFunction > const &function, double offset, double firstPulseWeight, double secondPulseWeight)
MANTID_CURVEFITTING_DLL Mantid::API::IFunction_sptr getDoublePulseMultiDomainFunction(std::shared_ptr< const API::MultiDomainFunction > const &function, double offset, double firstPulseWeight, double secondPulseWeight)
MANTID_CURVEFITTING_DLL Mantid::API::IFunction_sptr extractInnerFunction(std::shared_ptr< const Mantid::CurveFitting::Functions::Convolution > const &function)
std::unique_ptr< T > create(const P &parent, const IndexArg &indexArg, const HistArg &histArg)
This is the create() method that all the other create() methods call.
std::shared_ptr< IValidator > IValidator_sptr
A shared_ptr to an IValidator.
Definition: IValidator.h:26
std::string to_string(const wide_integer< Bits, Signed > &n)
@ InOut
Both an input & output workspace.
Definition: Property.h:55
@ Input
An input workspace.
Definition: Property.h:53
@ Output
An output workspace.
Definition: Property.h:54