Mantid
Loading...
Searching...
No Matches
LoadEventNexusIndexSetup.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 +
7#include <utility>
8
14#include "MantidIndexing/Extract.h"
15#include "MantidIndexing/Scatter.h"
16#include "MantidIndexing/SpectrumIndexSet.h"
17#include "MantidIndexing/SpectrumNumber.h"
18#include "MantidTypes/SpectrumDefinition.h"
19
20#include <algorithm>
21
22using namespace Mantid::API;
23using namespace Mantid::Indexing;
24
25namespace Mantid::DataHandling {
26
27namespace {
28void setupConsistentSpectrumNumbers(IndexInfo &filtered, const std::vector<detid_t> &detIDs) {
29 std::vector<Indexing::SpectrumNumber> spectrumNumbers;
30 // Temporary spectrum number in `filtered` was detector ID, now translate
31 // to spectrum number, starting at 1. Note that we use detIDs and not
32 // DetectorInfo for translation since we need to match the unfiltered
33 // spectrum numbers, which are based on skipping monitors (which would be
34 // included in DetectorInfo).
35 for (int32_t i = 0; i < static_cast<int32_t>(detIDs.size()); ++i) {
36 if (filtered.spectrumNumber(spectrumNumbers.size()) == detIDs[i])
37 spectrumNumbers.emplace_back(i + 1);
38 if (filtered.size() == spectrumNumbers.size())
39 break;
40 }
41 filtered.setSpectrumNumbers(std::move(spectrumNumbers));
42}
43} // namespace
44
46 const int32_t max, std::vector<int32_t> range,
47 Parallel::Communicator communicator)
48 : m_instrumentWorkspace(std::move(instrumentWorkspace)), m_min(min), m_max(max), m_range(std::move(range)),
49 m_communicator(std::move(communicator)) {}
50
51std::pair<int32_t, int32_t> LoadEventNexusIndexSetup::eventIDLimits() const { return {m_min, m_max}; }
52
54 // The default 1:1 will suffice but exclude the monitors as they are always in
55 // a separate workspace
56 auto detIDs = m_instrumentWorkspace->getInstrument()->getDetectorIDs(true);
57 const auto &detectorInfo = m_instrumentWorkspace->detectorInfo();
58 std::vector<SpectrumDefinition> specDefs;
59 specDefs.reserve(detIDs.size());
60 std::transform(detIDs.cbegin(), detIDs.cend(), std::back_inserter(specDefs),
61 [&detectorInfo](const auto detID) { return SpectrumDefinition(detectorInfo.indexOf(detID)); });
62 // We need to filter based on detector IDs, but use IndexInfo for filtering
63 // for a unified filtering mechanism. Thus we set detector IDs as (temporary)
64 // spectrum numbers.
65 IndexInfo indexInfo(std::vector<SpectrumNumber>(detIDs.begin(), detIDs.end()), Parallel::StorageMode::Cloned,
67 indexInfo.setSpectrumDefinitions(specDefs);
68
69 auto filtered = filterIndexInfo(indexInfo);
70
71 // Spectrum numbers are continuous and start at 1. If there is a filter,
72 // spectrum numbers are set up to be consistent with the unfiltered case.
73 if (filtered.size() == indexInfo.size()) {
74 filtered.setSpectrumNumbers(1, static_cast<int32_t>(filtered.size()));
75 } else {
76 setupConsistentSpectrumNumbers(filtered, detIDs);
77 }
78
79 return scatter(filtered);
80}
81
82IndexInfo LoadEventNexusIndexSetup::makeIndexInfo(const std::vector<std::string> &bankNames) {
83 const auto &componentInfo = m_instrumentWorkspace->componentInfo();
84 const auto &detectorInfo = m_instrumentWorkspace->detectorInfo();
85 std::vector<SpectrumDefinition> spectrumDefinitions;
86 // Temporary spectrum numbers setup up to be detector IDs, used for finding
87 // correct spectrum number to be consistent with unfiltered case.
88 std::vector<SpectrumNumber> spectrumNumbers;
89 const auto &instrument = m_instrumentWorkspace->getInstrument();
90 for (const auto &bankName : bankNames) {
91 const auto &bank = instrument->getComponentByName(bankName);
92 std::vector<size_t> dets;
93 if (bank) {
94 const auto bankIndex = componentInfo.indexOf(bank->getComponentID());
95 dets = componentInfo.detectorsInSubtree(bankIndex);
96 for (const auto detIndex : dets) {
97 spectrumDefinitions.emplace_back(detIndex);
98 spectrumNumbers.emplace_back(detectorInfo.detectorIDs()[detIndex]);
99 }
100 }
101 if (dets.empty())
102 throw std::runtime_error("Could not find the bank named '" + bankName +
103 "' as a component assembly in the instrument "
104 "tree; or it did not contain any detectors. Try "
105 "unchecking SingleBankPixelsOnly.");
106 }
107 Indexing::IndexInfo indexInfo(std::move(spectrumNumbers), Parallel::StorageMode::Cloned, m_communicator);
108 indexInfo.setSpectrumDefinitions(std::move(spectrumDefinitions));
109 setupConsistentSpectrumNumbers(indexInfo, instrument->getDetectorIDs(true));
110 // Filters are ignored when selecting bank names. Reset min/max to avoid
111 // unintended dropping of events in the loader.
112 m_min = EMPTY_INT();
113 m_max = EMPTY_INT();
114 return scatter(indexInfo);
115}
116
118 const std::pair<std::vector<int32_t>, std::vector<int32_t>> &spectrumDetectorMapping, const bool monitorsOnly) {
119 const auto &spec = spectrumDetectorMapping.first;
120 const auto &udet = spectrumDetectorMapping.second;
121
122 const std::vector<detid_t> monitors = m_instrumentWorkspace->getInstrument()->getMonitors();
123 const auto &detectorInfo = m_instrumentWorkspace->detectorInfo();
124 if (monitorsOnly) {
125 std::vector<Indexing::SpectrumNumber> spectrumNumbers;
126 std::vector<SpectrumDefinition> spectrumDefinitions;
127 // Find the det_ids in the udet array.
128 for (const auto id : monitors) {
129 // Find the index in the udet array
130 auto it = std::find(udet.begin(), udet.end(), id);
131 if (it != udet.end()) {
132 const specnum_t &specNo = spec[it - udet.begin()];
133 spectrumNumbers.emplace_back(specNo);
134 spectrumDefinitions.emplace_back(detectorInfo.indexOf(id));
135 }
136 }
137 Indexing::IndexInfo indexInfo(spectrumNumbers, Parallel::StorageMode::Cloned, m_communicator);
138 indexInfo.setSpectrumDefinitions(std::move(spectrumDefinitions));
139 return scatter(indexInfo);
140 } else {
141 SpectrumDetectorMapping mapping(spec, udet, monitors);
142 auto uniqueSpectra = mapping.getSpectrumNumbers();
143 std::vector<SpectrumDefinition> spectrumDefinitions;
144 for (const auto specNo : uniqueSpectra) {
145 spectrumDefinitions.emplace_back();
146 for (const auto detID : mapping.getDetectorIDsForSpectrumNo(specNo)) {
147 try {
148 spectrumDefinitions.back().add(detectorInfo.indexOf(detID));
149 } catch (std::out_of_range &) {
150 // Discarding detector IDs that do not exist in the instrument.
151 }
152 }
153 }
154 Indexing::IndexInfo indexInfo(std::vector<Indexing::SpectrumNumber>(uniqueSpectra.begin(), uniqueSpectra.end()),
155 Parallel::StorageMode::Cloned, m_communicator);
156 indexInfo.setSpectrumDefinitions(std::move(spectrumDefinitions));
157 return scatter(filterIndexInfo(indexInfo));
158 }
159}
160
166IndexInfo LoadEventNexusIndexSetup::filterIndexInfo(const IndexInfo &indexInfo) {
167 // Check if range [SpectrumMin, SpectrumMax] was supplied
168 if (m_min != EMPTY_INT() || m_max != EMPTY_INT()) {
169 if (m_max == EMPTY_INT())
170 m_max = static_cast<int32_t>(indexInfo.spectrumNumber(indexInfo.size() - 1));
171 if (m_min == EMPTY_INT())
172 m_min = static_cast<int32_t>(indexInfo.spectrumNumber(0));
173 // Avoid adding non-existing indices (can happen if instrument has gaps in
174 // its detector IDs). IndexInfo does the filtering for use.
175 const auto indices = indexInfo.makeIndexSet(static_cast<SpectrumNumber>(m_min), static_cast<SpectrumNumber>(m_max));
176 std::transform(indices.begin(), indices.end(), std::back_inserter(m_range),
177 [indexInfo](const auto index) { return static_cast<int32_t>(indexInfo.spectrumNumber(index)); });
178 }
179 // Check if SpectrumList was supplied (or filled via min/max above)
180 if (!m_range.empty()) {
181 std::sort(m_range.begin(), m_range.end());
182 const auto indices = indexInfo.makeIndexSet(std::vector<SpectrumNumber>(m_range.begin(), m_range.end()));
183 m_min = static_cast<int32_t>(indexInfo.spectrumNumber(*indices.begin()));
184 m_max = static_cast<int32_t>(indexInfo.spectrumNumber(*(indices.end() - 1)));
185 return extract(indexInfo, indices);
186 }
187 return indexInfo;
188}
189
190} // namespace Mantid::DataHandling
std::map< DeltaEMode::Type, std::string > index
Definition: DeltaEMode.cpp:19
specnum_t m_min
specnum_t m_max
A minimal class to hold the mapping between the spectrum number and its related detector ID numbers f...
const std::set< detid_t > & getDetectorIDsForSpectrumNo(const specnum_t spectrumNo) const
std::set< specnum_t > getSpectrumNumbers() const
LoadEventNexusIndexSetup(API::MatrixWorkspace_const_sptr instrumentWorkspace, const int32_t min, const int32_t max, std::vector< int32_t > range, Parallel::Communicator communicator=Parallel::Communicator())
const API::MatrixWorkspace_const_sptr m_instrumentWorkspace
Indexing::IndexInfo filterIndexInfo(const Indexing::IndexInfo &indexInfo)
Filter IndexInfo based on optional spectrum range/list provided.
std::pair< int32_t, int32_t > eventIDLimits() const
std::shared_ptr< const MatrixWorkspace > MatrixWorkspace_const_sptr
shared pointer to the matrix workspace base class (const version)
constexpr int EMPTY_INT() noexcept
Returns what we consider an "empty" integer within a property.
Definition: EmptyValues.h:25
int32_t specnum_t
Typedef for a spectrum Number.
Definition: IDTypes.h:16
STL namespace.