LoadMLZ.cpp 13.4 KB
Newer Older
Celine Durniak's avatar
Celine Durniak committed
1
#include "MantidDataHandling/LoadMLZ.h"
2
#include "MantidAPI/Axis.h"
3
#include "MantidGeometry/Instrument/DetectorInfo.h"
Celine Durniak's avatar
Celine Durniak committed
4
#include "MantidAPI/FileProperty.h"
5
#include "MantidAPI/MatrixWorkspace.h"
Celine Durniak's avatar
Celine Durniak committed
6
7
#include "MantidAPI/Progress.h"
#include "MantidAPI/RegisterFileLoader.h"
8
#include "MantidAPI/SpectrumInfo.h"
9
#include "MantidAPI/WorkspaceFactory.h"
10
#include "MantidDataHandling/LoadHelper.h"
Celine Durniak's avatar
Celine Durniak committed
11
12
#include "MantidGeometry/Instrument.h"
#include "MantidKernel/EmptyValues.h"
13
#include "MantidKernel/Exception.h"
Hahn, Steven's avatar
Hahn, Steven committed
14
#include "MantidKernel/OptionalBool.h"
Celine Durniak's avatar
Celine Durniak committed
15
16
17
18
#include "MantidKernel/UnitFactory.h"

#include <algorithm>
#include <cmath>
19
20
#include <limits>
#include <vector>
Celine Durniak's avatar
Celine Durniak committed
21

22
23
24
25
26
27
namespace Mantid {
namespace DataHandling {

using namespace Kernel;
using namespace API;
using namespace NeXus;
LamarMoore's avatar
LamarMoore committed
28
29
using HistogramData::BinEdges;
using HistogramData::Counts;
30
31
32
33
34
35

// Register the algorithm into the AlgorithmFactory
DECLARE_NEXUS_FILELOADER_ALGORITHM(LoadMLZ)

/** Constructor
 */
36
LoadMLZ::LoadMLZ()
37
38
39
40
41
    : m_numberOfTubes{0}, m_numberOfPixelsPerTube{0}, m_numberOfChannels{0},
      m_numberOfHistograms{0}, m_monitorElasticPeakPosition{0},
      m_wavelength{0.0}, m_channelWidth{0.0}, m_timeOfFlightDelay{0.0},
      m_monitorCounts{0}, m_chopper_speed{0.0}, m_chopper_ratio{0}, m_l1{0.0},
      m_l2{0.0}, m_t1{0.0}, m_supportedInstruments{"TOFTOF", "DNS"} {}
42
43
44
45
46
47
48
49

/// Algorithm's name for identification. @see Algorithm::name
const std::string LoadMLZ::name() const { return "LoadMLZ"; }

/// Algorithm's version for identification. @see Algorithm::version
int LoadMLZ::version() const { return 1; }

/// Algorithm's category for identification. @see Algorithm::category
50
const std::string LoadMLZ::category() const { return "DataHandling\\Nexus"; }
51
52
53
54

/** Initialize the algorithm's properties.
 */
void LoadMLZ::init() {
55
56
57
  const std::vector<std::string> exts{".nxs", ".hdf", ".hd5"};
  declareProperty(Kernel::make_unique<FileProperty>("Filename", "",
                                                    FileProperty::Load, exts),
58
59
                  "File path of the Data file to load");

60
61
62
  declareProperty(make_unique<WorkspaceProperty<>>("OutputWorkspace", "",
                                                   Direction::Output),
                  "The name to use for the output workspace");
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
}

/** Execute the algorithm.
 */
void LoadMLZ::exec() {
  // Retrieve filename
  std::string filenameData = getPropertyValue("Filename");

  // open the root node
  NeXus::NXRoot dataRoot(filenameData);
  NXEntry dataFirstEntry = dataRoot.openFirstEntry();

  loadInstrumentDetails(dataFirstEntry);
  loadTimeDetails(dataFirstEntry);

  initWorkSpace(dataFirstEntry);

  runLoadInstrument(); // just to get IDF contents
  initInstrumentSpecific();

  loadDataIntoTheWorkSpace(dataFirstEntry);
  loadRunDetails(dataFirstEntry);
  loadExperimentDetails(dataFirstEntry);
  maskDetectors(dataFirstEntry);

  // load the instrument from the IDF if it exists
  runLoadInstrument();
  setProperty("OutputWorkspace", m_localWorkspace);
}

/**
* Return the confidence with which this algorithm can load the file
* @param descriptor A descriptor for the file
* @returns An integer specifying the confidence level. 0 indicates it will not
* be used
*/
int LoadMLZ::confidence(Kernel::NexusDescriptor &descriptor) const {
  // fields existent only at the MLZ
  if (descriptor.pathExists("/Scan/wavelength") &&
      descriptor.pathExists("/Scan/title") &&
      descriptor.pathExists("/Scan/mode")) {
    return 80;
  } else {
    return 0;
Celine Durniak's avatar
Celine Durniak committed
107
  }
108
109
110
111
112
113
114
115
116
117
118
119
}

/**
 * Loads Masked detectors from the /Scan/instrument/Detector/pixel_mask
 */
void LoadMLZ::maskDetectors(NeXus::NXEntry &entry) {
  // path to the pixel_mask
  std::string pmpath = "instrument/detector/pixel_mask";

  NeXus::NXInt pmdata = entry.openNXInt(pmpath);
  // load the counts from the file into memory
  pmdata.load();
120
  g_log.debug() << "PMdata size: " << pmdata.size() << '\n';
121
122
123
  std::vector<int> masked_detectors(pmdata(), pmdata() + pmdata.size());

  g_log.debug() << "Number of masked detectors: " << masked_detectors.size()
124
                << '\n';
125

126
127
  auto &detInfo = m_localWorkspace->mutableDetectorInfo();
  std::vector<size_t> indicesToMask;
Hahn, Steven's avatar
Hahn, Steven committed
128
  for (auto masked_detector : masked_detectors) {
129
    g_log.debug() << "List of masked detectors: ";
130
    g_log.debug() << masked_detector;
131
    g_log.debug() << ", ";
132
133
134
135
136
137
    try {
      indicesToMask.push_back(detInfo.indexOf(masked_detector));
    } catch (std::out_of_range &) {
      g_log.warning() << "Invalid detector ID " << masked_detector
                      << ". Found while running LoadMLZ\n";
    }
Celine Durniak's avatar
Celine Durniak committed
138
  }
139
  g_log.debug() << '\n';
140

141
  for (const auto index : indicesToMask)
142
    detInfo.setMasked(index, true);
143
}
Celine Durniak's avatar
Celine Durniak committed
144

145
146
147
148
/**
 * Set the instrument name along with its path on the nexus file
*/
void LoadMLZ::loadInstrumentDetails(NeXus::NXEntry &firstEntry) {
Celine Durniak's avatar
Celine Durniak committed
149

150
  m_instrumentPath = m_mlzloader.findInstrumentNexusPath(firstEntry);
Celine Durniak's avatar
Celine Durniak committed
151

152
  if (m_instrumentPath.empty()) {
153
154
    throw std::runtime_error(
        "Cannot set the instrument name from the Nexus file!");
Celine Durniak's avatar
Celine Durniak committed
155
156
  }

157
158
  m_instrumentName = m_mlzloader.getStringFromNexusPath(
      firstEntry, m_instrumentPath + "/name");
Celine Durniak's avatar
Celine Durniak committed
159

160
161
162
163
164
  if (std::find(m_supportedInstruments.begin(), m_supportedInstruments.end(),
                m_instrumentName) == m_supportedInstruments.end()) {
    std::string message =
        "The instrument " + m_instrumentName + " is not valid for this loader!";
    throw std::runtime_error(message);
Celine Durniak's avatar
Celine Durniak committed
165
166
  }

167
  g_log.debug() << "Instrument name set to: " + m_instrumentName << '\n';
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
}

/**
 * Creates the workspace and initialises member variables with
 * the corresponding values
 *
 * @param entry :: The Nexus entry
 *
 */
void LoadMLZ::initWorkSpace(
    NeXus::NXEntry &entry) //, const std::vector<std::vector<int> >&monitors)
{
  // read in the data
  NXData dataGroup = entry.openNXData("data");
  NXInt data = dataGroup.openIntData();

  m_numberOfTubes = static_cast<size_t>(data.dim0());
  m_numberOfPixelsPerTube = static_cast<size_t>(data.dim1());
  m_numberOfChannels = static_cast<size_t>(data.dim2());
  m_numberOfHistograms = m_numberOfTubes * m_numberOfPixelsPerTube;

189
190
191
  g_log.debug() << "NumberOfTubes: " << m_numberOfTubes << '\n';
  g_log.debug() << "NumberOfPixelsPerTube: " << m_numberOfPixelsPerTube << '\n';
  g_log.debug() << "NumberOfChannels: " << m_numberOfChannels << '\n';
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206

  // Now create the output workspace
  m_localWorkspace = WorkspaceFactory::Instance().create(
      "Workspace2D", m_numberOfHistograms, m_numberOfChannels + 1,
      m_numberOfChannels);
  m_localWorkspace->getAxis(0)->unit() = UnitFactory::Instance().create("TOF");
  m_localWorkspace->setYUnitLabel("Counts");
}

/**
 * Function to do specific instrument stuff
 *
 */
void LoadMLZ::initInstrumentSpecific() {
  // Read data from IDF: distance source-sample and distance sample-detectors
207
208
  m_l1 = m_localWorkspace->spectrumInfo().l1();
  m_l2 = m_localWorkspace->spectrumInfo().l2(1);
209

210
  g_log.debug() << "L1: " << m_l1 << ", L2: " << m_l2 << '\n';
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
}

/**
 * Load the time details from the nexus file.
 * @param entry :: The Nexus entry
 */
void LoadMLZ::loadTimeDetails(NeXus::NXEntry &entry) {

  m_wavelength = entry.getFloat("wavelength");

  // Monitor can be monitor or Monitor
  std::string monitorName;
  if (entry.containsGroup("monitor"))
    monitorName = "monitor";
  else if (entry.containsGroup("Monitor"))
    monitorName = "Monitor";
  else {
    std::string message("Cannot find monitor/Monitor in the Nexus file!");
    g_log.error(message);
    throw std::runtime_error(message);
Celine Durniak's avatar
Celine Durniak committed
231
232
  }

233
  m_monitorCounts = entry.getInt(monitorName + "/integral");
Celine Durniak's avatar
Celine Durniak committed
234

235
  m_monitorElasticPeakPosition = entry.getInt(monitorName + "/elastic_peak");
Celine Durniak's avatar
Celine Durniak committed
236

237
238
239
  NXFloat time_of_flight_data =
      entry.openNXFloat(monitorName + "/time_of_flight");
  time_of_flight_data.load();
Celine Durniak's avatar
Celine Durniak committed
240

241
242
243
244
  // The entry "monitor/time_of_flight", has 3 fields:
  // channel width [microseconds], number of channels, Time of flight delay
  m_channelWidth = time_of_flight_data[0] * 50.e-3;
  m_timeOfFlightDelay = time_of_flight_data[2] * 50.e-3;
Celine Durniak's avatar
Celine Durniak committed
245

246
  g_log.debug("Nexus Data:");
247
248
249
  g_log.debug() << " MonitorCounts: " << m_monitorCounts << '\n';
  g_log.debug() << " ChannelWidth (microseconds): " << m_channelWidth << '\n';
  g_log.debug() << " Wavelength (angstroems): " << m_wavelength << '\n';
250
  g_log.debug() << " ElasticPeakPosition: " << m_monitorElasticPeakPosition
251
                << '\n';
252
  g_log.debug() << " TimeOfFlightDelay (microseconds): " << m_timeOfFlightDelay
253
                << '\n';
Celine Durniak's avatar
Celine Durniak committed
254

255
  m_chopper_speed = entry.getFloat("instrument/chopper/rotation_speed");
Celine Durniak's avatar
Celine Durniak committed
256

257
  m_chopper_ratio = entry.getInt("instrument/chopper/ratio");
Celine Durniak's avatar
Celine Durniak committed
258

259
260
  g_log.debug() << " ChopperSpeed: " << m_chopper_speed << '\n';
  g_log.debug() << " ChopperRatio: " << m_chopper_ratio << '\n';
261
}
Celine Durniak's avatar
Celine Durniak committed
262

263
264
265
266
267
268
269
270
/**
* Load information about the run.
* People from ISIS have this...
* TODO: They also have a lot of info in XML format!
*
* @param entry :: The Nexus entry
*/
void LoadMLZ::loadRunDetails(NXEntry &entry) {
Celine Durniak's avatar
Celine Durniak committed
271

272
  API::Run &runDetails = m_localWorkspace->mutableRun();
Celine Durniak's avatar
Celine Durniak committed
273

274
275
276
  std::string runNum = entry.getString("entry_identifier"); // run_number");
  std::string run_num = boost::lexical_cast<std::string>(runNum);
  runDetails.addProperty("run_number", run_num);
Celine Durniak's avatar
Celine Durniak committed
277

278
279
  std::string start_time = entry.getString("start_time");
  runDetails.addProperty("run_start", start_time);
Celine Durniak's avatar
Celine Durniak committed
280

281
282
  std::string end_time = entry.getString("end_time");
  runDetails.addProperty("run_end", end_time);
Celine Durniak's avatar
Celine Durniak committed
283

284
  runDetails.addProperty("wavelength", m_wavelength, "Angstrom", true);
Celine Durniak's avatar
Celine Durniak committed
285

286
  double ei = m_mlzloader.calculateEnergy(m_wavelength);
287
  runDetails.addProperty<double>("Ei", ei, "meV", true); // overwrite
Celine Durniak's avatar
Celine Durniak committed
288

289
290
  int duration = entry.getInt("duration");
  runDetails.addProperty("duration", duration, "Seconds", true);
Celine Durniak's avatar
Celine Durniak committed
291

292
293
  std::string mode = entry.getString("mode");
  runDetails.addProperty("mode", mode);
Celine Durniak's avatar
Celine Durniak committed
294

295
296
  std::string title = entry.getString("title");
  m_localWorkspace->setTitle(title);
Celine Durniak's avatar
Celine Durniak committed
297

298
299
300
  // Check if temperature is defined
  NXClass sample = entry.openNXGroup("sample");
  if (sample.containsDataSet("temperature")) {
301
302
    double temperature = entry.getFloat("sample/temperature");
    runDetails.addProperty("temperature", temperature, "K", true);
Celine Durniak's avatar
Celine Durniak committed
303
304
  }

305
306
307
308
  runDetails.addProperty("monitor_counts", m_monitorCounts);
  runDetails.addProperty("chopper_speed", m_chopper_speed);
  runDetails.addProperty("chopper_ratio", m_chopper_ratio);
  runDetails.addProperty("channel_width", m_channelWidth, "microseconds", true);
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328

  // Calculate number of full time channels - use to crop workspace - S. Busch's
  // method
  double full_channels = floor(30. * m_chopper_ratio / (m_chopper_speed)*1.e6 /
                               m_channelWidth); // channelWidth in microsec.
  runDetails.addProperty("full_channels", full_channels);

  // Proposal title
  std::string proposal_title = entry.getString("proposal");
  runDetails.addProperty("proposal_title", proposal_title);

  // proposal number
  std::string proposal_number = entry.getString("proposal_number");
  runDetails.addProperty("proposal_number", proposal_number);

  // users
  std::string user_name = entry.getString("user2/name");
  runDetails.addProperty("experiment_team", user_name);

  runDetails.addProperty("EPP", m_monitorElasticPeakPosition);
329
  runDetails.addProperty("TOF1", m_t1, "microseconds", true);
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
}

/**
 * Load data about the Experiment.
 *
 * TODO: This is very incomplete.
 *
 * @param entry :: The Nexus entry
 */
void LoadMLZ::loadExperimentDetails(NXEntry &entry) {
  // TODO: Do the rest
  // Pick out the geometry information

  std::string description =
      boost::lexical_cast<std::string>(entry.getFloat("sample/description"));

  m_localWorkspace->mutableSample().setName(description);
}

/**
 * Loads all the spectra into the workspace, including that from the monitor
 *
 * @param entry :: The Nexus entry
 */
void LoadMLZ::loadDataIntoTheWorkSpace(NeXus::NXEntry &entry) {
  // read in the data
  NXData dataGroup = entry.openNXData("data");
  NXInt data = dataGroup.openIntData();
  data.load();

Marina Ganeva's avatar
Marina Ganeva committed
360
  m_t1 = m_mlzloader.calculateTOF(m_l1, m_wavelength) * 1.0e+6;
361
  g_log.debug() << " t1 (microseconds): " << m_t1 << '\n';
362

363
364
365
  std::vector<double> detectorTofBins(m_numberOfChannels + 1);
  for (size_t i = 0; i < m_numberOfChannels + 1; ++i) {
    detectorTofBins[i] =
366
        m_channelWidth * static_cast<double>(static_cast<int>(i)) + m_t1 +
367
        m_channelWidth / 2;
Celine Durniak's avatar
Celine Durniak committed
368
369
  }

370
  // Assign calculated bins to first X axis
LamarMoore's avatar
LamarMoore committed
371
  BinEdges edges(std::move(detectorTofBins));
372

373
  Progress progress(this, 0.0, 1.0, m_numberOfTubes * m_numberOfPixelsPerTube);
374
375
376
377
378
379
  size_t spec = 0;
  for (size_t i = 0; i < m_numberOfTubes; ++i) {
    for (size_t j = 0; j < m_numberOfPixelsPerTube; ++j) {
      // Assign Y
      int *data_p = &data(static_cast<int>(i), static_cast<int>(j), 0);

LamarMoore's avatar
LamarMoore committed
380
381
      m_localWorkspace->setHistogram(
          spec, edges, Counts(data_p, data_p + m_numberOfChannels));
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397

      ++spec;
      progress.report();
    }
  }
}

/**
 * Run the Child Algorithm LoadInstrument.
 */
void LoadMLZ::runLoadInstrument() {
  IAlgorithm_sptr loadInst = createChildAlgorithm("LoadInstrument");

  // Now execute the Child Algorithm. Catch and log any error, but don't stop.
  try {
    loadInst->setPropertyValue("InstrumentName", m_instrumentName);
398
    g_log.debug() << "InstrumentName" << m_instrumentName << '\n';
399
    loadInst->setProperty<MatrixWorkspace_sptr>("Workspace", m_localWorkspace);
400
    loadInst->setProperty("RewriteSpectraMap",
401
                          Mantid::Kernel::OptionalBool(true));
402
403
404
405
406
    loadInst->execute();
  } catch (...) {
    g_log.information("Cannot load the instrument definition.");
  }
}
Celine Durniak's avatar
Celine Durniak committed
407
408
409

} // namespace DataHandling
} // namespace Mantid