Mantid
Loading...
Searching...
No Matches
MuonProcess.cpp
Go to the documentation of this file.
1// Mantid Repository : https://github.com/mantidproject/mantid
2//
3// Copyright © 2018 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 +
8
15
19
20// free functions
21namespace {
27bool isOne(int i) { return (i == 1); }
28} // namespace
29
31using namespace Kernel;
32using namespace API;
33using namespace DataObjects;
35
36// Register the algorithm into the AlgorithmFactory
37DECLARE_ALGORITHM(MuonProcess)
38
39//----------------------------------------------------------------------------------------------
41const std::string MuonProcess::name() const { return "MuonProcess"; }
42
44int MuonProcess::version() const { return 1; }
45
47const std::string MuonProcess::category() const { return "Workflow\\Muon"; }
48
49//----------------------------------------------------------------------------------------------
50
51//----------------------------------------------------------------------------------------------
52/*
53 * Initialize the algorithm's properties.
54 */
56 declareProperty(
57 std::make_unique<WorkspaceProperty<Workspace>>("InputWorkspace", "", Direction::Input, PropertyMode::Mandatory),
58 "Input workspace loaded from file (e.g. by LoadMuonNexus)");
59
60 std::vector<std::string> allowedModes{"CorrectAndGroup", "Analyse", "Combined"};
61 auto modeVal = std::make_shared<CompositeValidator>();
62 modeVal->add(std::make_shared<StringListValidator>(allowedModes));
63 modeVal->add(std::make_shared<MandatoryValidator<std::string>>());
64 declareProperty("Mode", "Combined", modeVal,
65 "Mode to run in. CorrectAndGroup applies dead time "
66 "correction and grouping; Analyse changes bin offset, "
67 "crops, rebins and calculates asymmetry; Combined does all "
68 "of the above.");
69
70 declareProperty(std::make_unique<ArrayProperty<int>>("SummedPeriodSet", Direction::Input),
71 "Comma-separated list of periods to be summed");
72
73 declareProperty(std::make_unique<ArrayProperty<int>>("SubtractedPeriodSet", Direction::Input),
74 "Comma-separated list of periods to be subtracted from the "
75 "SummedPeriodSet");
76
77 declareProperty("ApplyDeadTimeCorrection", false,
78 "Whether dead time correction should be applied to loaded workspace");
79 declareProperty(std::make_unique<WorkspaceProperty<TableWorkspace>>("DeadTimeTable", "", Direction::Input,
81 "Table with dead time information, e.g. from LoadMuonNexus."
82 "Must be specified if ApplyDeadTimeCorrection is set true.");
83 declareProperty(std::make_unique<WorkspaceProperty<TableWorkspace>>("DetectorGroupingTable", "", Direction::Input,
85 "Table with detector grouping information, e.g. from LoadMuonNexus.");
86
87 declareProperty("TimeZero", EMPTY_DBL(), "Value used for Time Zero correction");
88 declareProperty("LoadedTimeZero", EMPTY_DBL(), std::make_shared<MandatoryValidator<double>>(),
89 "Time Zero value loaded from file, e.g. from LoadMuonNexus.");
90 declareProperty(std::make_unique<ArrayProperty<double>>("RebinParams"),
91 "Params used for rebinning. If empty - rebinning is not done.");
92 declareProperty("Xmin", EMPTY_DBL(), "Minimal X value to include");
93 declareProperty("Xmax", EMPTY_DBL(), "Maximal X value to include");
94
95 std::vector<std::string> allowedTypes{"PairAsymmetry", "GroupAsymmetry", "GroupCounts"};
96 declareProperty("OutputType", "PairAsymmetry", std::make_shared<StringListValidator>(allowedTypes),
97 "What kind of workspace required for analysis.");
98
99 declareProperty("PairFirstIndex", EMPTY_INT(), "Workspace index of the first pair group");
100 declareProperty("PairSecondIndex", EMPTY_INT(), "Workspace index of the second pair group");
101 declareProperty("Alpha", 1.0, "Alpha value of the pair");
102
103 declareProperty("GroupIndex", EMPTY_INT(), "Workspace index of the group");
104
105 declareProperty(std::make_unique<WorkspaceProperty<Workspace>>("OutputWorkspace", "", Direction::Output),
106 "An output workspace.");
107
108 declareProperty("CropWorkspace", true,
109 "Determines if the input workspace "
110 "should be cropped at Xmax, Xmin is "
111 "still aplied.");
112
113 declareProperty("WorkspaceName", "", "The name of the input workspace");
114}
115
116//----------------------------------------------------------------------------------------------
121 Progress progress(this, 0.0, 1.0, 5);
122
123 // Supplied input workspace
124 Workspace_sptr inputWS = getProperty("InputWorkspace");
125
126 std::vector<int> summedPeriods = getProperty("SummedPeriodSet");
127 std::vector<int> subtractedPeriods = getProperty("SubtractedPeriodSet");
128 auto allPeriodsWS = std::make_shared<WorkspaceGroup>();
129 progress.report();
130
131 // Deal with single-period workspace
132 if (auto ws = std::dynamic_pointer_cast<MatrixWorkspace>(inputWS)) {
133 allPeriodsWS->addWorkspace(ws);
134 }
135 // Deal with multi-period workspace
136 else if (auto group = std::dynamic_pointer_cast<WorkspaceGroup>(inputWS)) {
137 allPeriodsWS = group;
138 }
139
140 progress.report();
141
142 // Check mode
143 const std::string mode = getProperty("Mode");
144
145 // Dead time correction and grouping
146 if (mode != "Analyse") {
147 bool applyDtc = getProperty("ApplyDeadTimeCorrection");
148 // Deal with dead time correction (if required)
149 if (applyDtc) {
150 TableWorkspace_sptr deadTimes = getProperty("DeadTimeTable");
151 allPeriodsWS = applyDTC(allPeriodsWS, deadTimes);
152 }
153 progress.report();
154 TableWorkspace_sptr grouping;
155 grouping = getProperty("DetectorGroupingTable");
156 progress.report();
157 // Deal with grouping
158 allPeriodsWS = groupWorkspaces(allPeriodsWS, grouping);
159 } else {
160 progress.report();
161 progress.report();
162 }
163
164 // If not analysing, the present WS will be the output
165 Workspace_sptr outWS = allPeriodsWS;
166 if (mode != "CorrectAndGroup") {
167 // Correct bin values
168 double loadedTimeZero = getProperty("LoadedTimeZero");
169 allPeriodsWS = correctWorkspaces(allPeriodsWS, loadedTimeZero);
170
171 // Perform appropriate calculation
172 std::string outputType = getProperty("OutputType");
173 int groupIndex = getProperty("GroupIndex");
174 std::unique_ptr<IMuonAsymmetryCalculator> asymCalc;
175 if (outputType == "GroupCounts") {
176 asymCalc =
177 std::make_unique<MuonGroupCountsCalculator>(allPeriodsWS, summedPeriods, subtractedPeriods, groupIndex);
178 } else if (outputType == "GroupAsymmetry") {
179 asymCalc = std::make_unique<MuonGroupAsymmetryCalculator>(allPeriodsWS, summedPeriods, subtractedPeriods,
180 groupIndex, getProperty("Xmin"), getProperty("Xmax"),
181 getProperty("WorkspaceName"));
182 } else if (outputType == "PairAsymmetry") {
183 int first = getProperty("PairFirstIndex");
184 int second = getProperty("PairSecondIndex");
185 double alpha = getProperty("Alpha");
186 asymCalc = std::make_unique<MuonPairAsymmetryCalculator>(allPeriodsWS, summedPeriods, subtractedPeriods, first,
187 second, alpha);
188 }
189 progress.report();
190 outWS = asymCalc->calculate();
191 }
192
193 setProperty("OutputWorkspace", outWS);
194}
195
204 const TableWorkspace_sptr &grouping) {
205 WorkspaceGroup_sptr outWS = std::make_shared<WorkspaceGroup>();
206 for (int i = 0; i < wsGroup->getNumberOfEntries(); i++) {
207 auto ws = std::dynamic_pointer_cast<MatrixWorkspace>(wsGroup->getItem(i));
208 if (ws) {
210 auto group = createChildAlgorithm("MuonGroupDetectors");
211 group->setProperty("InputWorkspace", ws);
212 group->setProperty("DetectorGroupingTable", grouping);
213 group->execute();
214 result = group->getProperty("OutputWorkspace");
215 outWS->addWorkspace(result);
216 }
217 }
218 return outWS;
219}
220
228 WorkspaceGroup_sptr outWS = std::make_shared<WorkspaceGroup>();
229 for (int i = 0; i < wsGroup->getNumberOfEntries(); i++) {
230 auto ws = std::dynamic_pointer_cast<MatrixWorkspace>(wsGroup->getItem(i));
231 if (ws) {
233 auto dtc = createChildAlgorithm("ApplyDeadTimeCorr");
234 dtc->setProperty("InputWorkspace", ws);
235 dtc->setProperty("DeadTimeTable", dt);
236 dtc->execute();
237 result = dtc->getProperty("OutputWorkspace");
238 if (result) {
239 outWS->addWorkspace(result);
240 } else {
241 throw std::runtime_error("ApplyDeadTimeCorr failed to apply dead time "
242 "correction in MuonProcess");
243 }
244 }
245 }
246 return outWS;
247}
248
258 WorkspaceGroup_sptr outWS = std::make_shared<WorkspaceGroup>();
259 for (int i = 0; i < wsGroup->getNumberOfEntries(); i++) {
260 auto ws = std::dynamic_pointer_cast<MatrixWorkspace>(wsGroup->getItem(i));
261 if (ws) {
263 result = correctWorkspace(ws, loadedTimeZero);
264 outWS->addWorkspace(result);
265 }
266 }
267 return outWS;
268}
269
278 // Offset workspace, if need to
279 double timeZero = getProperty("TimeZero");
280 if (timeZero != EMPTY_DBL()) {
281 double offset = loadedTimeZero - timeZero;
282
283 auto changeOffset = createChildAlgorithm("ChangeBinOffset");
284 changeOffset->setProperty("InputWorkspace", ws);
285 changeOffset->setProperty("Offset", offset);
286 changeOffset->execute();
287
288 ws = changeOffset->getProperty("OutputWorkspace");
289 }
290 // Crop workspace, if need to
291 double Xmin = getProperty("Xmin");
292 double Xmax = getProperty("Xmax");
293 if (Xmin != EMPTY_DBL() || Xmax != EMPTY_DBL()) {
294 auto crop = createChildAlgorithm("CropWorkspace");
295 crop->setProperty("InputWorkspace", ws);
296
297 if (Xmin != EMPTY_DBL())
298 crop->setProperty("Xmin", Xmin);
299 bool toCrop = getProperty("CropWorkspace");
300 if (toCrop) {
301 if (Xmax != EMPTY_DBL())
302 crop->setProperty("Xmax", Xmax);
303 }
304 crop->execute();
305
306 ws = crop->getProperty("OutputWorkspace");
307 }
308
309 // Rebin workspace if need to
310 std::vector<double> rebinParams = getProperty("RebinParams");
311 if (!rebinParams.empty()) {
312 auto rebin = createChildAlgorithm("Rebin");
313 rebin->setProperty("InputWorkspace", ws);
314 rebin->setProperty("Params", rebinParams);
315 rebin->setProperty("FullBinsOnly", true);
316 rebin->execute();
317
318 ws = rebin->getProperty("OutputWorkspace");
319 }
320
321 return ws;
322}
323
334std::map<std::string, std::string> MuonProcess::validateInputs() {
335 std::map<std::string, std::string> errors;
336
337 // Supplied input workspace and sets of periods
338 const std::string propInputWS("InputWorkspace"), propSummedPeriodSet("SummedPeriodSet"),
339 propSubtractedPeriodSet("SubtractedPeriodSet");
340 Workspace_sptr inputWS = getProperty(propInputWS);
341 std::vector<int> summedPeriods = getProperty(propSummedPeriodSet);
342 std::vector<int> subtractedPeriods = getProperty(propSubtractedPeriodSet);
343
344 // If single-period data, test the sets of periods specified
345 if (auto ws = std::dynamic_pointer_cast<MatrixWorkspace>(inputWS)) {
346 if (std::find_if_not(summedPeriods.begin(), summedPeriods.end(), isOne) != summedPeriods.end()) {
347 errors[propSummedPeriodSet] = "Single period data but set of periods to "
348 "sum contains invalid values.";
349 }
350 if (!subtractedPeriods.empty()) {
351 errors[propSubtractedPeriodSet] = "Single period data but second set of periods specified";
352 }
353 } else {
354 // If not a MatrixWorkspace, must be a multi-period WorkspaceGroup
355 auto group = std::dynamic_pointer_cast<WorkspaceGroup>(inputWS);
356 if (group == nullptr) {
357 errors[propInputWS] = "Input workspace is of invalid type";
358 } else {
359 auto numPeriods = group->getNumberOfEntries();
360 if (numPeriods < 1) {
361 errors[propInputWS] = "Input workspace contains no periods";
362 }
363 // check summed period numbers
364 std::vector<int> invalidPeriods;
365 std::copy_if(summedPeriods.cbegin(), summedPeriods.cend(), std::back_inserter(invalidPeriods),
366 [numPeriods](auto period) { return period < 1 || period > numPeriods; });
367 if (!invalidPeriods.empty()) {
368 errors[propSummedPeriodSet] = buildErrorString(invalidPeriods);
369 invalidPeriods.clear();
370 }
371 // check subtracted period numbers
372 std::copy_if(subtractedPeriods.cbegin(), subtractedPeriods.cend(), std::back_inserter(invalidPeriods),
373 [numPeriods](auto period) { return period < 1 || period > numPeriods; });
374 if (!invalidPeriods.empty()) {
375 errors[propSubtractedPeriodSet] = buildErrorString(invalidPeriods);
376 }
377 }
378 }
379
380 // Some parameters can be mandatory or not depending on the mode
381 const std::string propMode("Mode"), propApplyDTC("ApplyDeadTimeCorrection"), propDeadTime("DeadTimeTable"),
382 propDetGroup("DetectorGroupingTable");
383 const std::string mode = getProperty(propMode);
384 // If analysis will take place, SummedPeriodSet is mandatory
385 if (mode != "CorrectAndGroup") {
386 if (summedPeriods.empty()) {
387 errors[propSummedPeriodSet] = "Cannot analyse: list of periods to sum was empty";
388 }
389 }
390 // If correcting/grouping will take place, DetectorGroupingTable is mandatory
391 if (mode != "Analyse") {
392 TableWorkspace_sptr grouping = getProperty(propDetGroup);
393 if (grouping == nullptr) {
394 errors[propDetGroup] = "No detector grouping table supplied";
395 }
396 }
397 // If dead time correction is to be applied, must supply dead times
398 bool applyDtc = getProperty(propApplyDTC);
399 if (applyDtc) {
400 TableWorkspace_sptr deadTimes = getProperty(propDeadTime);
401 if (!deadTimes) {
402 errors[propDeadTime] = "Cannot apply dead time correction as no dead times were supplied";
403 }
404 }
405
406 return errors;
407}
408
414std::string MuonProcess::buildErrorString(const std::vector<int> &invalidPeriods) const {
415 std::stringstream message;
416 message << "Invalid periods specified: ";
417 for (auto it = invalidPeriods.begin(); it != invalidPeriods.end(); it++) {
418 message << *it;
419 if (it != invalidPeriods.end() - 1) {
420 message << ", ";
421 }
422 }
423 return message.str();
424}
425
426} // namespace Mantid::WorkflowAlgorithms
#define DECLARE_ALGORITHM(classname)
Definition: Algorithm.h:576
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) override
Create a Child Algorithm.
Kernel::IPropertyManager::TypedValue getProperty(const std::string &name) const override
Get the property held by this object.
Helper class for reporting progress from algorithms.
Definition: Progress.h:25
A property class for workspaces.
Support for a property that holds an array of values.
Definition: ArrayProperty.h:28
Validator to check that a property is not left empty.
void report()
Increments the loop counter by 1, then sends the progress notification on behalf of its algorithm.
Definition: ProgressBase.h:51
MuonProcess : Processes and analyses Muon workspace.
Definition: MuonProcess.h:18
API::MatrixWorkspace_sptr correctWorkspace(API::MatrixWorkspace_sptr ws, double loadedTimeZero)
Applies offset, crops and rebin the workspace according to specified params.
std::map< std::string, std::string > validateInputs() override
Perform validation of inputs to the algorithm.
API::WorkspaceGroup_sptr applyDTC(const API::WorkspaceGroup_sptr &wsGroup, const DataObjects::TableWorkspace_sptr &dt)
Applies dead time correction to the workspace group.
std::string buildErrorString(const std::vector< int > &invalidPeriods) const
Builds an error message from a list of invalid periods.
API::WorkspaceGroup_sptr correctWorkspaces(const API::WorkspaceGroup_sptr &wsGroup, double loadedTimeZero)
Applies offset, crops and rebins all workspaces in the group.
const std::string category() const override
Algorithm's category for identification.
Definition: MuonProcess.cpp:47
int version() const override
Algorithm's version for identification.
Definition: MuonProcess.cpp:44
void exec() override
Execute the algorithm.
API::WorkspaceGroup_sptr groupWorkspaces(const API::WorkspaceGroup_sptr &wsGroup, const DataObjects::TableWorkspace_sptr &grouping)
Groups specified workspace group according to specified DetectorGroupingTable.
std::shared_ptr< WorkspaceGroup > WorkspaceGroup_sptr
shared pointer to Mantid::API::WorkspaceGroup
std::shared_ptr< Workspace > Workspace_sptr
shared pointer to Mantid::API::Workspace
Definition: Workspace_fwd.h:20
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
std::shared_ptr< TableWorkspace > TableWorkspace_sptr
shared pointer to Mantid::DataObjects::TableWorkspace
constexpr int EMPTY_INT() noexcept
Returns what we consider an "empty" integer within a property.
Definition: EmptyValues.h:25
constexpr double EMPTY_DBL() noexcept
Returns what we consider an "empty" double within a property.
Definition: EmptyValues.h:43
STL namespace.
@ Input
An input workspace.
Definition: Property.h:53
@ Output
An output workspace.
Definition: Property.h:54