Mantid
Loading...
Searching...
No Matches
IntegrateEllipsoidsTwoStep.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
10#include "MantidAPI/Run.h"
11#include "MantidAPI/Sample.h"
22
27
28#include <boost/math/special_functions/round.hpp>
29#include <cmath>
30#include <string>
31#include <tuple>
32
33using namespace Mantid::API;
34using namespace Mantid::DataObjects;
35using namespace Mantid::Kernel;
36
37namespace Mantid::MDAlgorithms {
38
39// Register the algorithm into the AlgorithmFactory
40DECLARE_ALGORITHM(IntegrateEllipsoidsTwoStep)
41
42//---------------------------------------------------------------------
44const std::string IntegrateEllipsoidsTwoStep::name() const { return "IntegrateEllipsoidsTwoStep"; }
45
47int IntegrateEllipsoidsTwoStep::version() const { return 1; }
48
50const std::string IntegrateEllipsoidsTwoStep::category() const { return "Crystal\\Integration"; }
51
53 auto ws_valid = std::make_shared<CompositeValidator>();
54 ws_valid->add<InstrumentValidator>();
55
56 auto mustBePositive = std::make_shared<BoundedValidator<double>>();
57 mustBePositive->setLower(0.0);
58
60 std::make_unique<WorkspaceProperty<MatrixWorkspace>>("InputWorkspace", "", Direction::Input, ws_valid),
61 "An input MatrixWorkspace with time-of-flight units along "
62 "X-axis and defined instrument with defined sample");
63
64 declareProperty(std::make_unique<WorkspaceProperty<PeaksWorkspace>>("PeaksWorkspace", "", Direction::InOut),
65 "Workspace with peaks to be integrated");
66
67 declareProperty("RegionRadius", .35, mustBePositive,
68 "Only events at most this distance from a peak will be "
69 "considered when integrating");
70
71 declareProperty("SpecifySize", false, "If true, use the following for the major axis sizes, else use 3-sigma");
72
73 declareProperty("PeakSize", .18, mustBePositive, "Half-length of major axis for peak ellipsoid");
74
75 declareProperty("BackgroundInnerSize", .18, mustBePositive,
76 "Half-length of major axis for inner ellipsoidal surface of "
77 "background region");
78
79 declareProperty("BackgroundOuterSize", .23, mustBePositive,
80 "Half-length of major axis for outer ellipsoidal surface of "
81 "background region");
82
83 declareProperty("IntegrateInHKL", false, "If true, integrate in HKL space not Q space.");
84 declareProperty("IntegrateIfOnEdge", true,
85 "Set to false to not integrate if peak radius is off edge of detector."
86 "Background will be scaled if background radius is off edge.");
87
88 declareProperty("AdaptiveQBackground", false,
89 "Default is false. If true, "
90 "BackgroundOuterRadius + AdaptiveQMultiplier * **|Q|** and "
91 "BackgroundInnerRadius + AdaptiveQMultiplier * **|Q|**");
92
93 declareProperty("AdaptiveQMultiplier", 0.0,
94 "PeakRadius + AdaptiveQMultiplier * **|Q|** "
95 "so each peak has a "
96 "different integration radius. Q includes the 2*pi factor.");
97
98 declareProperty("WeakPeakThreshold", 1.0, mustBePositive, "Intensity threshold use to classify a peak as weak.");
99
100 declareProperty("UseOnePercentBackgroundCorrection", true,
101 "If this options is enabled, then the the top 1% of the "
102 "background will be removed"
103 "before the background subtraction.");
104
105 declareProperty(std::make_unique<WorkspaceProperty<PeaksWorkspace>>("OutputWorkspace", "", Direction::Output),
106 "The output PeaksWorkspace will be a copy of the input PeaksWorkspace "
107 "with the peaks' integrated intensities.");
108}
109
111 PeaksWorkspace_sptr input_peak_ws = getProperty("PeaksWorkspace");
112 MatrixWorkspace_sptr input_ws = getProperty("InputWorkspace");
113 EventWorkspace_sptr eventWS = std::dynamic_pointer_cast<EventWorkspace>(input_ws);
114
115 Workspace2D_sptr histoWS = std::dynamic_pointer_cast<Workspace2D>(input_ws);
116 if (!eventWS && !histoWS) {
117 throw std::runtime_error("IntegrateEllipsoids needs either a "
118 "EventWorkspace or Workspace2D as input.");
119 }
120
121 const double weakPeakThreshold = getProperty("WeakPeakThreshold");
122
123 // validation of inputs
124 if (!input_peak_ws) {
125 throw std::runtime_error("Could not read the Peaks Workspace");
126 }
127
128 if (!input_ws) {
129 throw std::runtime_error("Could not read the Input Workspace");
130 }
131
132 PeaksWorkspace_sptr peak_ws = getProperty("OutputWorkspace");
133 if (peak_ws != input_peak_ws) {
134 peak_ws = input_peak_ws->clone();
135 }
136
137 Progress prog(this, 0.5, 1.0, input_ws->getNumberHistograms());
138
139 std::vector<Peak> &peaks = peak_ws->getPeaks();
140 size_t n_peaks = peak_ws->getNumberPeaks();
141 size_t indexed_count = 0;
142 std::vector<V3D> peak_q_list;
143 std::vector<V3D> hkl_vectors;
144 for (size_t i = 0; i < n_peaks; i++) // Note: we skip un-indexed peaks
145 {
146 V3D hkl(peaks[i].getH(), peaks[i].getK(), peaks[i].getL());
147 if (Geometry::IndexingUtils::ValidIndex(hkl, 1.0)) // use tolerance == 1 to
148 // just check for (0,0,0)
149 {
150 peak_q_list.emplace_back(peaks[i].getQLabFrame());
151 V3D miller_ind(static_cast<double>(boost::math::iround<double>(hkl[0])),
152 static_cast<double>(boost::math::iround<double>(hkl[1])),
153 static_cast<double>(boost::math::iround<double>(hkl[2])));
154 hkl_vectors.emplace_back(miller_ind);
155 indexed_count++;
156 }
157 }
158
159 if (indexed_count < 3) {
160 throw std::runtime_error("At least three linearly independent indexed peaks are needed.");
161 }
162 // Get UB using indexed peaks and
163 // lab-Q vectors
164 Matrix<double> UB(3, 3, false);
165 Geometry::IndexingUtils::Optimize_UB(UB, hkl_vectors, peak_q_list);
166 Matrix<double> UBinv(UB);
167 UBinv.Invert();
168 UBinv *= (1.0 / (2.0 * M_PI));
169
170 std::vector<std::pair<std::pair<double, double>, V3D>> qList;
171 for (size_t i = 0; i < n_peaks; i++) {
172 qList.emplace_back(std::pair<double, double>(1.0, 1.0), V3D(peaks[i].getQLabFrame()));
173 }
174
175 const bool integrateEdge = getProperty("IntegrateIfOnEdge");
176 if (!integrateEdge) {
177 // This only fails in the unit tests which say that MaskBTP is not
178 // registered
179 try {
180 runMaskDetectors(input_peak_ws, "Tube", "edges");
181 runMaskDetectors(input_peak_ws, "Pixel", "edges");
182 } catch (...) {
183 g_log.error("Can't execute MaskBTP algorithm for this instrument to set "
184 "edge for IntegrateIfOnEdge option");
185 }
186 calculateE1(input_peak_ws->detectorInfo()); // fill E1Vec for use in detectorQ
187 }
188
189 const bool integrateInHKL = getProperty("IntegrateInHKL");
190 bool useOnePercentBackgroundCorrection = getProperty("UseOnePercentBackgroundCorrection");
191 Integrate3DEvents integrator(qList, UBinv, getProperty("RegionRadius"), useOnePercentBackgroundCorrection);
192
193 if (eventWS) {
194 // process as EventWorkspace
195 qListFromEventWS(integrator, prog, eventWS, UBinv, integrateInHKL);
196 } else {
197 // process as Workspace2D
198 qListFromHistoWS(integrator, prog, histoWS, UBinv, integrateInHKL);
199 }
200
201 std::vector<std::pair<int, V3D>> weakPeaks, strongPeaks;
202
203 // Compute signal to noise ratio for all peaks
204 for (int index = 0; static_cast<size_t>(index) < qList.size(); ++index) {
205 const auto &item = qList[index];
206 const auto center = item.second;
208 auto sig2noise = integrator.estimateSignalToNoiseRatio(params, center);
209
210 auto &peak = peak_ws->getPeak(index);
211 peak.setIntensity(0);
212 peak.setSigmaIntensity(0);
213
214 const auto result = std::make_pair(index, center);
215 if (sig2noise < weakPeakThreshold) {
216 g_log.notice() << "Peak " << peak.getHKL() << " with Q = " << center << " is a weak peak with signal to noise "
217 << sig2noise << "\n";
218 weakPeaks.emplace_back(result);
219 } else {
220 g_log.notice() << "Peak " << peak.getHKL() << " with Q = " << center << " is a strong peak with signal to noise "
221 << sig2noise << "\n";
222 strongPeaks.emplace_back(result);
223 }
224 }
225
226 std::vector<std::pair<std::shared_ptr<const Geometry::PeakShape>, std::tuple<double, double, double>>> shapeLibrary;
227
228 // Integrate strong peaks
229 for (const auto &item : strongPeaks) {
230 const auto index = item.first;
231 const auto &q = item.second;
232 double inti, sigi;
233
235 const auto result = integrator.integrateStrongPeak(params, q, inti, sigi);
236 shapeLibrary.emplace_back(result);
237
238 auto &peak = peak_ws->getPeak(index);
239 peak.setIntensity(inti);
240 peak.setSigmaIntensity(sigi);
241 peak.setPeakShape(std::get<0>(result));
242 }
243
244 std::vector<Eigen::Vector3d> points;
245 std::transform(strongPeaks.begin(), strongPeaks.end(), std::back_inserter(points),
246 [&](const std::pair<int, V3D> &item) {
247 const auto q = item.second;
248 return Eigen::Vector3d(q[0], q[1], q[2]);
249 });
250
251 if (points.empty())
252 throw std::runtime_error("Cannot integrate peaks when all peaks are below "
253 "the signal to noise ratio.");
254
255 NearestNeighbours<3> kdTree(points);
256
257 // Integrate weak peaks
258 for (const auto &item : weakPeaks) {
259 double inti, sigi;
260 const auto index = item.first;
261 const auto &q = item.second;
262
263 const auto result = kdTree.findNearest(Eigen::Vector3d(q[0], q[1], q[2]));
264 const auto strongIndex = static_cast<int>(std::get<1>(result[0]));
265
266 auto &peak = peak_ws->getPeak(index);
267 auto &strongPeak = peak_ws->getPeak(strongIndex);
268
269 g_log.notice() << "Integrating weak peak " << peak.getHKL() << " using strong peak " << strongPeak.getHKL() << "\n";
270
271 const auto libShape = shapeLibrary[static_cast<int>(strongIndex)];
272 const auto shape = std::dynamic_pointer_cast<const PeakShapeEllipsoid>(libShape.first);
273 const auto frac = std::get<0>(libShape.second);
274
275 g_log.notice() << "Weak peak will be adjusted by " << frac << "\n";
276 IntegrationParameters params = makeIntegrationParameters(strongPeak.getQLabFrame());
277 const auto weakShape = integrator.integrateWeakPeak(params, shape, libShape.second, q, inti, sigi);
278
279 peak.setIntensity(inti);
280 peak.setSigmaIntensity(sigi);
281 peak.setPeakShape(weakShape);
282 }
283
284 // This flag is used by the PeaksWorkspace to evaluate whether it has been
285 // integrated.
286 peak_ws->mutableRun().addProperty("PeaksIntegrated", 1, true);
287 setProperty("OutputWorkspace", peak_ws);
288}
289
292 params.peakRadius = getProperty("PeakSize");
293 params.backgroundInnerRadius = getProperty("BackgroundInnerSize");
294 params.backgroundOuterRadius = getProperty("BackgroundOuterSize");
295 params.regionRadius = getProperty("RegionRadius");
296 params.specifySize = getProperty("SpecifySize");
297 params.E1Vectors = E1Vec;
298
299 const bool adaptiveQBackground = getProperty("AdaptiveQBackground");
300 const double adaptiveQMultiplier = getProperty("AdaptiveQMultiplier");
301 const double adaptiveQBackgroundMultiplier = (adaptiveQBackground) ? adaptiveQMultiplier : 0.0;
302
303 // modulus of Q
304 const double lenQpeak = peak_q.norm();
305 // change params to support adaptive Q
306 params.peakRadius = adaptiveQMultiplier * lenQpeak + params.peakRadius;
307 params.backgroundInnerRadius = adaptiveQBackgroundMultiplier * lenQpeak + params.backgroundInnerRadius;
308 params.backgroundOuterRadius = adaptiveQBackgroundMultiplier * lenQpeak + params.backgroundOuterRadius;
309 return params;
310}
311
313 EventWorkspace_sptr &wksp, DblMatrix const &UBinv, bool hkl_integ) {
314 // loop through the eventlists
315
316 const std::string ELASTIC("Elastic");
318 const std::string Q3D("Q3D");
319 const std::size_t DIMS(3);
320
321 MDWSDescription m_targWSDescr;
322 m_targWSDescr.setMinMax(std::vector<double>(3, -2000.), std::vector<double>(3, 2000.));
323 m_targWSDescr.buildFromMatrixWS(wksp, Q3D, ELASTIC);
324 m_targWSDescr.setLorentsCorr(false);
325
326 // generate the detectors table
327 Mantid::API::Algorithm_sptr childAlg = createChildAlgorithm("PreprocessDetectorsToMD", 0.,
328 .5); // HACK. soft dependency on non-dependent package.
329 childAlg->setProperty("InputWorkspace", wksp);
330 childAlg->executeAsChildAlg();
331
332 DataObjects::TableWorkspace_sptr table = childAlg->getProperty("OutputWorkspace");
333 if (!table)
334 throw(std::runtime_error("Can not retrieve results of \"PreprocessDetectorsToMD\""));
335
336 m_targWSDescr.m_PreprDetTable = table;
337
338 auto numSpectra = static_cast<int>(wksp->getNumberHistograms());
340 for (int i = 0; i < numSpectra; ++i) {
342
343 // units conversion helper
344 UnitsConversionHelper unitConverter;
345 unitConverter.initialize(m_targWSDescr, "Momentum");
346
347 // initialize the MD coordinates conversion class
348 MDTransfQ3D qConverter;
349 qConverter.initialize(m_targWSDescr);
350
351 std::vector<double> buffer(DIMS);
352 // get a reference to the event list
353 EventList &events = wksp->getSpectrum(i);
354
356 events.compressEvents(1e-5, &events);
357
358 // check to see if the event list is empty
359 if (events.empty()) {
360 prog.report();
361 continue; // nothing to do
362 }
363
364 // update which pixel is being converted
365 std::vector<Mantid::coord_t> locCoord(DIMS, 0.);
366 unitConverter.updateConversion(i);
367 qConverter.calcYDepCoordinates(locCoord, i);
368
369 // loop over the events
370 double signal(1.); // ignorable garbage
371 double errorSq(1.); // ignorable garbage
372 const std::vector<WeightedEventNoTime> &raw_events = events.getWeightedEventsNoTime();
373 std::vector<std::pair<std::pair<double, double>, V3D>> qList;
374 for (const auto &raw_event : raw_events) {
375 double val = unitConverter.convertUnits(raw_event.tof());
376 qConverter.calcMatrixCoord(val, locCoord, signal, errorSq);
377 for (size_t dim = 0; dim < DIMS; ++dim) {
378 buffer[dim] = locCoord[dim];
379 }
380 V3D qVec(buffer[0], buffer[1], buffer[2]);
381 if (hkl_integ)
382 qVec = UBinv * qVec;
383 qList.emplace_back(std::pair<double, double>(raw_event.m_weight, raw_event.m_errorSquared), qVec);
384 } // end of loop over events in list
385 PARALLEL_CRITICAL(addEvents) { integrator.addEvents(qList, hkl_integ); }
386
387 prog.report();
389 } // end of loop over spectra
391}
392
403 DblMatrix const &UBinv, bool hkl_integ) {
404
405 // loop through the eventlists
406 const std::string ELASTIC("Elastic");
408 const std::string Q3D("Q3D");
409 const std::size_t DIMS(3);
410
411 MDWSDescription m_targWSDescr;
412 m_targWSDescr.setMinMax(std::vector<double>(3, -2000.), std::vector<double>(3, 2000.));
413 m_targWSDescr.buildFromMatrixWS(wksp, Q3D, ELASTIC);
414 m_targWSDescr.setLorentsCorr(false);
415
416 // generate the detectors table
417 Mantid::API::Algorithm_sptr childAlg = createChildAlgorithm("PreprocessDetectorsToMD", 0.,
418 .5); // HACK. soft dependency on non-dependent package.
419 childAlg->setProperty("InputWorkspace", wksp);
420 childAlg->executeAsChildAlg();
421
422 DataObjects::TableWorkspace_sptr table = childAlg->getProperty("OutputWorkspace");
423 if (!table)
424 throw(std::runtime_error("Can not retrieve results of \"PreprocessDetectorsToMD\""));
425 else
426 m_targWSDescr.m_PreprDetTable = table;
427
428 auto numSpectra = static_cast<int>(wksp->getNumberHistograms());
430 for (int i = 0; i < numSpectra; ++i) {
432
433 // units conversion helper
434 UnitsConversionHelper unitConverter;
435 unitConverter.initialize(m_targWSDescr, "Momentum");
436
437 // initialize the MD coordinates conversion class
438 MDTransfQ3D qConverter;
439 qConverter.initialize(m_targWSDescr);
440
441 // get tof and y values
442 const auto &xVals = wksp->points(i);
443 const auto &yVals = wksp->y(i);
444 const auto &eVals = wksp->e(i);
445
446 // update which pixel is being converted
447 std::vector<Mantid::coord_t> locCoord(DIMS, 0.);
448 unitConverter.updateConversion(i);
449 qConverter.calcYDepCoordinates(locCoord, i);
450
451 // loop over the events
452 double signal(1.); // ignorable garbage
453 double errorSq(1.); // ignorable garbage
454
455 std::vector<std::pair<std::pair<double, double>, V3D>> qList;
456
457 for (size_t j = 0; j < yVals.size(); ++j) {
458 const double &yVal = yVals[j];
459 const double &esqVal = eVals[j] * eVals[j]; // error squared (variance)
460 if (yVal > 0) // TODO, is this condition right?
461 {
462 double val = unitConverter.convertUnits(xVals[j]);
463 qConverter.calcMatrixCoord(val, locCoord, signal, errorSq);
464 V3D qVec(locCoord[0], locCoord[1], locCoord[2]);
465 if (hkl_integ)
466 qVec = UBinv * qVec;
467
468 if (std::isnan(qVec[0]) || std::isnan(qVec[1]) || std::isnan(qVec[2]))
469 continue;
470 // Account for counts in histograms by increasing the qList with the
471 // same q-point
472 qList.emplace_back(std::pair<double, double>(yVal, esqVal), qVec);
473 }
474 }
475 PARALLEL_CRITICAL(addHisto) { integrator.addEvents(qList, hkl_integ); }
476 prog.report();
478 } // end of loop over spectra
480}
481
482/*
483 * Define edges for each instrument by masking. For CORELLI, tubes 1 and 16, and
484 *pixels 0 and 255.
485 * Get Q in the lab frame for every peak, call it C
486 * For every point on the edge, the trajectory in reciprocal space is a straight
487 *line, going through O=V3D(0,0,0).
488 * Calculate a point at a fixed momentum, say k=1. Q in the lab frame
489 *E=V3D(-k*sin(tt)*cos(ph),-k*sin(tt)*sin(ph),k-k*cos(ph)).
490 * Normalize E to 1: E=E*(1./E.norm())
491 *
492 * @param inst: instrument
493 */
495 for (size_t i = 0; i < detectorInfo.size(); ++i) {
496 if (detectorInfo.isMonitor(i))
497 continue; // skip monitor
498 if (!detectorInfo.isMasked(i))
499 continue; // edge is masked so don't check if not masked
500 const auto &det = detectorInfo.detector(i);
501 double tt1 = det.getTwoTheta(V3D(0, 0, 0), V3D(0, 0, 1)); // two theta
502 double ph1 = det.getPhi(); // phi
503 V3D E1 = V3D(-std::sin(tt1) * std::cos(ph1), -std::sin(tt1) * std::sin(ph1),
504 1. - std::cos(tt1)); // end of trajectory
505 E1 = E1 * (1. / E1.norm()); // normalize
506 E1Vec.emplace_back(E1);
507 }
508}
509
511 const std::string &property, const std::string &values) {
512 auto alg = createChildAlgorithm("MaskBTP");
513 alg->setProperty<Workspace_sptr>("Workspace", peakWS);
514 alg->setProperty(property, values);
515 if (!alg->execute())
516 throw std::runtime_error("MaskDetectors Child Algorithm has not executed successfully");
517}
518} // namespace Mantid::MDAlgorithms
#define DECLARE_ALGORITHM(classname)
Definition: Algorithm.h:576
std::map< DeltaEMode::Type, std::string > index
Definition: DeltaEMode.cpp:19
#define PARALLEL_START_INTERRUPT_REGION
Begins a block to skip processing is the algorithm has been interupted Note the end of the block if n...
Definition: MultiThreaded.h:94
#define PARALLEL_CRITICAL(name)
#define PARALLEL_END_INTERRUPT_REGION
Ends a block to skip processing is the algorithm has been interupted Note the start of the block if n...
#define PARALLEL_FOR_IF(condition)
Empty definitions - to enable set your complier to enable openMP.
#define PARALLEL_CHECK_INTERRUPT_REGION
Adds a check after a Parallel region to see if it was interupted.
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
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
Kernel::Logger & g_log
Definition: Algorithm.h:451
A validator which checks that a workspace has a valid instrument.
Helper class for reporting progress from algorithms.
Definition: Progress.h:25
A property class for workspaces.
A class for holding :
Definition: EventList.h:56
void compressEvents(double tolerance, EventList *destination)
Compress the event list by grouping events with the same TOF (within a given tolerance).
Definition: EventList.cpp:1676
std::vector< WeightedEventNoTime > & getWeightedEventsNoTime()
Return the list of WeightedEvent contained.
Definition: EventList.cpp:823
void switchTo(Mantid::API::EventType newType) override
Switch the EventList to use the given EventType (TOF, WEIGHTED, or WEIGHTED_NOTIME)
Definition: EventList.cpp:649
bool empty() const
Much like stl containers, returns true if there is nothing in the event list.
Definition: EventList.cpp:1158
Geometry::DetectorInfo is an intermediate step towards a DetectorInfo that is part of Instrument-2....
Definition: DetectorInfo.h:49
bool isMasked(const size_t index) const
Returns true if the detector is masked.
const Geometry::IDetector & detector(const size_t index) const
Return a const reference to the detector with given index.
size_t size() const
Returns the size of the DetectorInfo, i.e., the number of detectors in the instrument.
bool isMonitor(const size_t index) const
Returns true if the detector is a monitor.
virtual double getTwoTheta(const Kernel::V3D &observer, const Kernel::V3D &axis) const =0
Gives the angle of this detector object with respect to an axis.
static double Optimize_UB(Kernel::DblMatrix &UB, const std::vector< Kernel::V3D > &hkl_vectors, const std::vector< Kernel::V3D > &q_vectors, std::vector< double > &sigabc)
Find the UB matrix that most nearly maps hkl to qxyz for 3 or more peaks.
static bool ValidIndex(const Kernel::V3D &hkl, double tolerance)
Check is hkl is within tolerance of integer (h,k,l) non-zero values.
IPropertyManager * setProperty(const std::string &name, const T &value)
Templated method to set the value of a PropertyWithValue.
void notice(const std::string &msg)
Logs at notice level.
Definition: Logger.cpp:95
void error(const std::string &msg)
Logs at error level.
Definition: Logger.cpp:77
Numerical Matrix class.
Definition: Matrix.h:42
T Invert()
LU inversion routine.
Definition: Matrix.cpp:924
NearestNeighbourResults findNearest(const VectorType &pos, const size_t k=1, const double error=0.0)
Find the k nearest neighbours to a given point.
void report()
Increments the loop counter by 1, then sends the progress notification on behalf of its algorithm.
Definition: ProgressBase.h:51
Class for 3D vectors.
Definition: V3D.h:34
double norm() const noexcept
Definition: V3D.h:263
This is a low-level class to construct a map with lists of events near each peak Q-vector,...
std::pair< std::shared_ptr< const Mantid::Geometry::PeakShape >, std::tuple< double, double, double > > integrateStrongPeak(const IntegrationParameters &params, const Kernel::V3D &peak_q, double &inti, double &sigi)
Find the net integrated intensity of a peak, using ellipsoidal volumes.
std::shared_ptr< const Geometry::PeakShape > integrateWeakPeak(const IntegrationParameters &params, Mantid::DataObjects::PeakShapeEllipsoid_const_sptr shape, const std::tuple< double, double, double > &libPeak, const Mantid::Kernel::V3D &peak_q, double &inti, double &sigi)
double estimateSignalToNoiseRatio(const IntegrationParameters &params, const Mantid::Kernel::V3D &center, bool forceSpherical=false, double sphericityTol=0.02)
void addEvents(std::vector< std::pair< std::pair< double, double >, Mantid::Kernel::V3D > > const &event_qs, bool hkl_integ)
Add event Q's to lists of events near peaks.
IntegrateEllipsoidsTwoStep provides a two pass peak integration algorithm.
void qListFromHistoWS(Integrate3DEvents &integrator, API::Progress &prog, DataObjects::Workspace2D_sptr &wksp, const Kernel::DblMatrix &UBinv, bool hkl_integ)
qListFromHistoWS creates qlist from input workspaces of type Workspace2D
void exec() override
Virtual method - must be overridden by concrete algorithm.
void init() override
Virtual method - must be overridden by concrete algorithm.
IntegrationParameters makeIntegrationParameters(const Kernel::V3D &peak_q) const
void qListFromEventWS(Integrate3DEvents &integrator, API::Progress &prog, DataObjects::EventWorkspace_sptr &wksp, const Kernel::DblMatrix &UBinv, bool hkl_integ)
int version() const override
Get the version of this algorithm.
const std::string category() const override
Get the category of this algorithm.
void runMaskDetectors(const Mantid::DataObjects::PeaksWorkspace_sptr &peakWS, const std::string &property, const std::string &values)
void calculateE1(const Geometry::DetectorInfo &detectorInfo)
Calculate if this Q is on a detector.
std::vector< Kernel::V3D > E1Vec
save all detector pixels
Class responsible for conversion of input workspace data into proper number of output dimensions for ...
Definition: MDTransfQ3D.h:31
bool calcMatrixCoord(const double &deltaEOrK0, std::vector< coord_t > &Coord, double &s, double &err) const override
Calculates 3D transformation of the variable coordinates and (if applicable) signal and error dependi...
Definition: MDTransfQ3D.cpp:47
void initialize(const MDWSDescription &ConvParams) override
function initalizes all variables necessary for converting workspace variables into MD variables in M...
bool calcYDepCoordinates(std::vector< coord_t > &Coord, size_t i) override
Method updates the value of preprocessed detector coordinates in Q-space, used by other functions.
helper class describes the properties of target MD workspace, which should be obtained as the result ...
void setMinMax(const std::vector< double > &minVal, const std::vector< double > &maxVal)
function sets up min-max values to the dimensions, described by the class
void buildFromMatrixWS(const API::MatrixWorkspace_sptr &pWS, const std::string &QMode, const std::string &dEMode, const std::vector< std::string > &dimPropertyNames=std::vector< std::string >())
method builds MD Event ws description from a matrix workspace and the transformations,...
DataObjects::TableWorkspace_const_sptr m_PreprDetTable
void setLorentsCorr(bool On=false)
do we need to perform Lorentz corrections
void initialize(const MDWSDescription &targetWSDescr, const std::string &unitsTo, bool forceViaTOF=false)
Initialize unit conversion helper This method is interface to internal initialize method,...
std::complex< double > MANTID_API_DLL E1(std::complex< double > z)
Integral for Gamma.
std::shared_ptr< Workspace > Workspace_sptr
shared pointer to Mantid::API::Workspace
Definition: Workspace_fwd.h:20
std::shared_ptr< Algorithm > Algorithm_sptr
Typedef for a shared pointer to an Algorithm.
Definition: Algorithm.h:61
@ WEIGHTED_NOTIME
Definition: IEventList.h:18
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
std::shared_ptr< Workspace2D > Workspace2D_sptr
shared pointer to Mantid::DataObjects::Workspace2D
std::shared_ptr< TableWorkspace > TableWorkspace_sptr
shared pointer to Mantid::DataObjects::TableWorkspace
std::shared_ptr< PeaksWorkspace > PeaksWorkspace_sptr
Typedef for a shared pointer to a peaks workspace.
std::shared_ptr< EventWorkspace > EventWorkspace_sptr
shared pointer to the EventWorkspace class
std::enable_if< std::is_pointer< Arg >::value, bool >::type threadSafe(Arg workspace)
Thread-safety check Checks the workspace to ensure it is suitable for multithreaded access.
Definition: MultiThreaded.h:22
const std::string Q3D("Q3D")
Only convert to Q-vector.
const std::size_t DIMS(3)
Q-vector is always three dimensional.
const std::string ELASTIC("Elastic")
This only works for diffraction.
STL namespace.
@ 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