-
Notifications
You must be signed in to change notification settings - Fork 122
/
CreateLogTimeCorrection.cpp
158 lines (128 loc) · 6.07 KB
/
CreateLogTimeCorrection.cpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
// Mantid Repository : https://github.com/mantidproject/mantid
//
// Copyright © 2018 ISIS Rutherford Appleton Laboratory UKRI,
// NScD Oak Ridge National Laboratory, European Spallation Source,
// Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
// SPDX - License - Identifier: GPL - 3.0 +
#include "MantidAlgorithms/CreateLogTimeCorrection.h"
#include "MantidAPI/FileProperty.h"
#include "MantidAPI/InstrumentValidator.h"
#include "MantidAPI/MatrixWorkspace.h"
#include "MantidAPI/Run.h"
#include "MantidAPI/TableRow.h"
#include "MantidAPI/WorkspaceProperty.h"
#include "MantidGeometry/Instrument/DetectorInfo.h"
#include <fstream>
using namespace Mantid;
using namespace Mantid::API;
using namespace Mantid::DataObjects;
using namespace Mantid::Geometry;
using namespace Mantid::Kernel;
using namespace std;
namespace Mantid::Algorithms {
DECLARE_ALGORITHM(CreateLogTimeCorrection)
//----------------------------------------------------------------------------------------------
/** Declare properties
*/
void CreateLogTimeCorrection::init() {
declareProperty(std::make_unique<WorkspaceProperty<MatrixWorkspace>>("InputWorkspace", "", Direction::Input,
std::make_shared<InstrumentValidator>()),
"Name of the input workspace to generate log correct from.");
declareProperty(std::make_unique<WorkspaceProperty<TableWorkspace>>("OutputWorkspace", "", Direction::Output),
"Name of the output workspace containing the corrections.");
declareProperty(std::make_unique<FileProperty>("OutputFilename", "", FileProperty::OptionalSave),
"Name of the output time correction file.");
}
//----------------------------------------------------------------------------------------------
/** Main execution body
*/
void CreateLogTimeCorrection::exec() {
// 1. Process input
MatrixWorkspace_sptr dataWS = getProperty("InputWorkspace");
// Check whether the output workspace name is same as input
string outwsname = getPropertyValue("OutputWorkspace");
if (outwsname == dataWS->getName()) {
stringstream errmsg;
errmsg << "It is not allowed to use the same name by both input matrix "
"workspace and output table workspace.";
g_log.error(errmsg.str());
throw runtime_error(errmsg.str());
}
const auto &detectorInfo = dataWS->detectorInfo();
// 2. Log some information
logGeometryInformation(detectorInfo);
// 3. Calculate log time correction
auto corrections = calculateCorrections(detectorInfo);
// 4. Output
TableWorkspace_sptr outWS = generateCorrectionTable(detectorInfo, corrections);
setProperty("OutputWorkspace", outWS);
string filename = getProperty("OutputFilename");
g_log.information() << "Output file name is " << filename << ".\n";
if (!filename.empty()) {
writeCorrectionToFile(filename, detectorInfo, corrections);
}
}
//----------------------------------------------------------------------------------------------
/** Get instrument geometry setup including L2 for each detector and L1
*/
void CreateLogTimeCorrection::logGeometryInformation(const Geometry::DetectorInfo &detectorInfo) const {
g_log.information() << "Sample position = " << detectorInfo.samplePosition() << "; "
<< "Source position = " << detectorInfo.sourcePosition() << ", L1 = " << detectorInfo.l1() << "; "
<< "Number of detector/pixels = " << detectorInfo.size() << ".\n";
}
//----------------------------------------------------------------------------------------------
/** Calculate the log time correction for each pixel, i.e., correcton from event
* time at detector
* to time at sample
*/
std::vector<double> CreateLogTimeCorrection::calculateCorrections(const Geometry::DetectorInfo &detectorInfo) const {
std::vector<double> corrections(detectorInfo.size());
const double l1 = detectorInfo.l1();
for (size_t detectorIndex = 0; detectorIndex < detectorInfo.size(); ++detectorIndex) {
double corrfactor = l1 / (l1 + detectorInfo.l2(detectorIndex));
corrections[detectorIndex] = corrfactor;
}
return corrections;
}
//----------------------------------------------------------------------------------------------
/** Write L2 map and correction map to a TableWorkspace
*/
TableWorkspace_sptr CreateLogTimeCorrection::generateCorrectionTable(const Geometry::DetectorInfo &detectorInfo,
const std::vector<double> &corrections) const {
auto tablews = std::make_shared<TableWorkspace>();
tablews->addColumn("int", "DetectorID");
tablews->addColumn("double", "Correction");
tablews->addColumn("double", "L2");
const auto &detectorIds = detectorInfo.detectorIDs();
for (size_t detectorIndex = 0; detectorIndex < detectorInfo.size(); ++detectorIndex) {
if (!detectorInfo.isMonitor(detectorIndex)) {
const detid_t detid = detectorIds[detectorIndex];
const double correction = corrections[detectorIndex];
const double l2 = detectorInfo.l2(detectorIndex);
TableRow newrow = tablews->appendRow();
newrow << detid << correction << l2;
}
}
return tablews;
}
//----------------------------------------------------------------------------------------------
/** Write correction map to a text file
*/
void CreateLogTimeCorrection::writeCorrectionToFile(const string &filename, const Geometry::DetectorInfo &detectorInfo,
const std::vector<double> &corrections) const {
ofstream ofile;
ofile.open(filename.c_str());
if (ofile.is_open()) {
const auto &detectorIds = detectorInfo.detectorIDs();
for (size_t detectorIndex = 0; detectorIndex < corrections.size(); ++detectorIndex) {
if (!detectorInfo.isMonitor(detectorIndex)) {
ofile << detectorIds[detectorIndex] << "\t" << setw(20) << setprecision(5) << corrections[detectorIndex]
<< "\n";
}
}
ofile.close();
} else {
g_log.error() << "Unable to open file " << filename << " to write!\n";
}
}
} // namespace Mantid::Algorithms