LoadILLSANS.cpp 39.3 KB
Newer Older
1
2
3
// Mantid Repository : https://github.com/mantidproject/mantid
//
// Copyright © 2018 ISIS Rutherford Appleton Laboratory UKRI,
4
5
//   NScD Oak Ridge National Laboratory, European Spallation Source,
//   Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
6
// SPDX - License - Identifier: GPL - 3.0 +
7
#include "MantidDataHandling/LoadILLSANS.h"
8
#include "MantidAPI/Axis.h"
9
#include "MantidAPI/FileProperty.h"
10
#include "MantidAPI/MatrixWorkspace.h"
11
#include "MantidAPI/Progress.h"
12
#include "MantidAPI/RegisterFileLoader.h"
13
#include "MantidAPI/SpectrumInfo.h"
14
#include "MantidAPI/WorkspaceFactory.h"
15
#include "MantidDataHandling/LoadHelper.h"
16
#include "MantidGeometry/IDetector.h"
Hahn, Steven's avatar
Hahn, Steven committed
17
#include "MantidGeometry/Instrument.h"
18
#include "MantidGeometry/Instrument/RectangularDetector.h"
Hahn, Steven's avatar
Hahn, Steven committed
19
#include "MantidHistogramData/LinearGenerator.h"
20
#include "MantidKernel/BoundedValidator.h"
21
#include "MantidKernel/ConfigService.h"
22
#include "MantidKernel/EmptyValues.h"
Hahn, Steven's avatar
Hahn, Steven committed
23
#include "MantidKernel/OptionalBool.h"
24
#include "MantidKernel/PhysicalConstants.h"
25
#include "MantidKernel/UnitFactory.h"
26
#include "MantidKernel/VectorHelper.h"
27

LamarMoore's avatar
LamarMoore committed
28
#include <Poco/Path.h>
29

30
namespace Mantid::DataHandling {
31

Gagik Vardanyan's avatar
Gagik Vardanyan committed
32
33
34
35
namespace {
static constexpr size_t N_MONITORS = 2;
}

36
37
38
39
using namespace Kernel;
using namespace API;
using namespace NeXus;

40
DECLARE_NEXUS_FILELOADER_ALGORITHM(LoadILLSANS)
41
42
43
44

//----------------------------------------------------------------------------------------------
/** Constructor
 */
45
LoadILLSANS::LoadILLSANS()
Samuel Jones's avatar
Samuel Jones committed
46
47
    : m_supportedInstruments{"D11", "D22", "D33", "D16"}, m_defaultBinning{0, 0}, m_resMode("nominal"), m_isTOF(false),
      m_sourcePos(0.) {}
48
49
50

//----------------------------------------------------------------------------------------------
/// Algorithm's name for identification. @see Algorithm::name
51
const std::string LoadILLSANS::name() const { return "LoadILLSANS"; }
52
53

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

/// Algorithm's category for identification. @see Algorithm::category
Samuel Jones's avatar
Samuel Jones committed
57
const std::string LoadILLSANS::category() const { return "DataHandling\\Nexus;ILL\\SANS"; }
58

59
60
/// Algorithm's summary. @see Algorithm::summery
const std::string LoadILLSANS::summary() const {
61
  return "Loads ILL nexus files for SANS instruments D11, D16, D22, D33.";
62
63
}

64
65
66
67
68
//----------------------------------------------------------------------------------------------

/**
 * Return the confidence with with this algorithm can load the file
 * @param descriptor A descriptor for the file
69
70
 * @returns An integer specifying the confidence level. 0 indicates it will not
 * be used
71
 */
72
73
int LoadILLSANS::confidence(Kernel::NexusDescriptor &descriptor) const {
  // fields existent only at the ILL for SANS machines
74
  if (descriptor.pathExists("/entry0/mode") &&
Samuel Jones's avatar
Samuel Jones committed
75
76
77
78
      ((descriptor.pathExists("/entry0/reactor_power") && descriptor.pathExists("/entry0/instrument_name")) ||
       (descriptor.pathExists("/entry0/instrument/name") && descriptor.pathExists("/entry0/acquisition_mode") &&
        !descriptor.pathExists("/entry0/instrument/Detector")))) // serves to remove the TOF
                                                                 // instruments
79
  {
80
    return 80;
81
82
83
  } else {
    return 0;
  }
84
85
86
87
88
89
}

//----------------------------------------------------------------------------------------------
/** Initialize the algorithm's properties.
 */
void LoadILLSANS::init() {
Samuel Jones's avatar
Samuel Jones committed
90
  declareProperty(std::make_unique<FileProperty>("Filename", "", FileProperty::Load, ".nxs"),
Sam Jenkins's avatar
Sam Jenkins committed
91
                  "Name of the nexus file to load");
Samuel Jones's avatar
Samuel Jones committed
92
  declareProperty(std::make_unique<WorkspaceProperty<>>("OutputWorkspace", "", Direction::Output),
93
                  "The name to use for the output workspace");
94
95
  auto mustBePositive = std::make_shared<BoundedValidator<double>>();
  mustBePositive->setLower(0);
Samuel Jones's avatar
Samuel Jones committed
96
97
98
  declareProperty("Wavelength", 0.0, mustBePositive,
                  "The wavelength of the experiment, in angstroms. Used only for D16. Will "
                  "override the nexus' value if there is one.");
99
100
101
102
103
104
}

//----------------------------------------------------------------------------------------------
/** Execute the algorithm.
 */
void LoadILLSANS::exec() {
105
  const std::string filename = getPropertyValue("Filename");
106
  m_isD16Omega = false;
107
108
  NXRoot root(filename);
  NXEntry firstEntry = root.openFirstEntry();
Samuel Jones's avatar
Samuel Jones committed
109
  const std::string instrumentPath = m_loadHelper.findInstrumentNexusPath(firstEntry);
110
  setInstrumentName(firstEntry, instrumentPath);
111
  Progress progress(this, 0.0, 1.0, 4);
112
  progress.report("Initializing the workspace for " + m_instrumentName);
113
114
  if (m_instrumentName == "D33") {
    initWorkSpaceD33(firstEntry, instrumentPath);
115
    progress.report("Loading the instrument " + m_instrumentName);
116
    runLoadInstrument();
Samuel Jones's avatar
Samuel Jones committed
117
    const DetectorPosition detPos = getDetectorPositionD33(firstEntry, instrumentPath);
118
    progress.report("Moving detectors");
119
    moveDetectorsD33(detPos);
Gagik Vardanyan's avatar
Gagik Vardanyan committed
120
    if (m_isTOF) {
Gagik Vardanyan's avatar
Gagik Vardanyan committed
121
122
      adjustTOF();
      moveSource();
Gagik Vardanyan's avatar
Gagik Vardanyan committed
123
    }
124

125
126
127
128
129
  } else if (m_instrumentName == "D16") {
    initWorkSpace(firstEntry, instrumentPath);
    progress.report("Loading the instrument " + m_instrumentName);
    runLoadInstrument();

Samuel Jones's avatar
Samuel Jones committed
130
    double distance = firstEntry.getFloat(instrumentPath + "/Det/value") / 1000; // mm to metre
131
    const double angle = firstEntry.getFloat(instrumentPath + "/Gamma/value");
132
    placeD16(-angle, distance, "detector");
133

Mathieu Tillet's avatar
Mathieu Tillet committed
134
135
136
137
  } else if (m_instrumentName == "D11B") {
    initWorkSpaceD11B(firstEntry, instrumentPath);
    progress.report("Loading the instrument " + m_instrumentName);
    runLoadInstrument();
138
139
140

    // we move the parent "detector" component, but since it is at (0,0,0), we
    // need to find the distance it has to move and move it to this position
Samuel Jones's avatar
Samuel Jones committed
141
    double finalDistance = firstEntry.getFloat(instrumentPath + "/Detector 1/det_calc");
142
143
    V3D pos = getComponentPosition("detector_center");
    double currentDistance = pos.Z();
144

145
146
147
    moveDetectorDistance(finalDistance - currentDistance, "detector");
    API::Run &runDetails = m_localWorkspace->mutableRun();
    runDetails.addProperty<double>("L2", finalDistance, true);
Mathieu Tillet's avatar
Mathieu Tillet committed
148
149
150

  } else if (m_instrumentName == "D22B") {
    initWorkSpaceD22B(firstEntry, instrumentPath);
151

152
153
154
155
156
    const std::string backIndex = m_localWorkspace->getInstrument()->getStringParameter("back_detector_index")[0];
    const std::string frontIndex = m_localWorkspace->getInstrument()->getStringParameter("front_detector_index")[0];

    // first we move the central (back) detector
    double distance = firstEntry.getFloat(instrumentPath + "/Detector " + backIndex + "/det" + backIndex + "_calc");
157
    moveDetectorDistance(distance, "detector_back");
158
159
    API::Run &runDetails = m_localWorkspace->mutableRun();
    runDetails.addProperty<double>("L2", distance, true);
160

161
    double offset = firstEntry.getFloat(instrumentPath + "/Detector " + backIndex + "/dtr" + backIndex + "_actual");
162
    moveDetectorHorizontal(-offset / 1000, "detector_back"); // mm to meter
163

164
165
    // then the front (right) one
    distance = firstEntry.getFloat(instrumentPath + "/Detector " + frontIndex + "/det" + frontIndex + "_calc");
166
    moveDetectorDistance(distance, "detector_front");
Mathieu Tillet's avatar
Mathieu Tillet committed
167

168
    // mm to meter
169
    offset = firstEntry.getFloat(instrumentPath + "/Detector " + frontIndex + "/dtr" + frontIndex + "_actual");
170
    moveDetectorHorizontal(-offset / 1000, "detector_front"); // mm to meter
171
    double angle = firstEntry.getFloat(instrumentPath + "/Detector " + frontIndex + "/dan" + frontIndex + "_actual");
172
    rotateInstrument(-angle, "detector_front");
173

174
  } else {
175
    // D11 and D22
176
    initWorkSpace(firstEntry, instrumentPath);
177
    progress.report("Loading the instrument " + m_instrumentName);
178
    runLoadInstrument();
Samuel Jones's avatar
Samuel Jones committed
179
    double distance = m_loadHelper.getDoubleFromNexusPath(firstEntry, instrumentPath + "/detector/det_calc");
180
    progress.report("Moving detectors");
181
    moveDetectorDistance(distance, "detector");
182
183
    API::Run &runDetails = m_localWorkspace->mutableRun();
    runDetails.addProperty<double>("L2", distance, true);
184
    if (m_instrumentName == "D22") {
Samuel Jones's avatar
Samuel Jones committed
185
      double offset = m_loadHelper.getDoubleFromNexusPath(firstEntry, instrumentPath + "/detector/dtr_actual");
186
      moveDetectorHorizontal(-offset / 1000, "detector"); // mm to meter
187
    }
188
  }
189

190
  progress.report("Setting sample logs");
191
  setFinalProperties(filename);
192
  setProperty("OutputWorkspace", m_localWorkspace);
193
}
194
195
196

/**
 * Set member variable with the instrument name
197
198
 * @param firstEntry: already opened first entry in nexus
 * @param instrumentNamePath : the path inside nexus where the instrument name is written
199
 */
Samuel Jones's avatar
Samuel Jones committed
200
void LoadILLSANS::setInstrumentName(const NeXus::NXEntry &firstEntry, const std::string &instrumentNamePath) {
201
  if (instrumentNamePath.empty()) {
Samuel Jones's avatar
Samuel Jones committed
202
    throw std::runtime_error("Cannot set the instrument name from the Nexus file!");
203
  }
Samuel Jones's avatar
Samuel Jones committed
204
205
  m_instrumentName = m_loadHelper.getStringFromNexusPath(firstEntry, instrumentNamePath + "/name");
  const auto inst = std::find(m_supportedInstruments.begin(), m_supportedInstruments.end(), m_instrumentName);
Mathieu Tillet's avatar
Mathieu Tillet committed
206

Samuel Jones's avatar
Samuel Jones committed
207
  if ((m_instrumentName == "D11" || m_instrumentName == "D22") && firstEntry.containsGroup("data1")) {
Mathieu Tillet's avatar
Mathieu Tillet committed
208
    m_instrumentName += "B";
Mathieu Tillet's avatar
Mathieu Tillet committed
209
210
  }

211
  if (inst == m_supportedInstruments.end()) {
Samuel Jones's avatar
Samuel Jones committed
212
213
    throw std::runtime_error("Instrument " + m_instrumentName +
                             " is not supported. Only D11, D16, D22 and D33 are supported");
214
  }
215
  g_log.debug() << "Instrument name set to: " + m_instrumentName << '\n';
216
217
218
219
}

/**
 * Get detector panel distances from the nexus file
220
221
 * @param firstEntry : already opened first entry in nexus
 * @param instrumentNamePath : the path inside nexus where the instrument name is written
222
223
 * @return a structure with the positions
 */
Samuel Jones's avatar
Samuel Jones committed
224
225
LoadILLSANS::DetectorPosition LoadILLSANS::getDetectorPositionD33(const NeXus::NXEntry &firstEntry,
                                                                  const std::string &instrumentNamePath) {
226
227
  std::string detectorPath(instrumentNamePath + "/detector");
  DetectorPosition pos;
Samuel Jones's avatar
Samuel Jones committed
228
229
230
231
232
233
234
235
  pos.distanceSampleRear = m_loadHelper.getDoubleFromNexusPath(firstEntry, detectorPath + "/det2_calc");
  pos.distanceSampleBottomTop = m_loadHelper.getDoubleFromNexusPath(firstEntry, detectorPath + "/det1_calc");
  pos.distanceSampleRightLeft = pos.distanceSampleBottomTop + m_loadHelper.getDoubleFromNexusPath(
                                                                  firstEntry, detectorPath + "/det1_panel_separation");
  pos.shiftLeft = m_loadHelper.getDoubleFromNexusPath(firstEntry, detectorPath + "/OxL_actual") * 1e-3;
  pos.shiftRight = m_loadHelper.getDoubleFromNexusPath(firstEntry, detectorPath + "/OxR_actual") * 1e-3;
  pos.shiftUp = m_loadHelper.getDoubleFromNexusPath(firstEntry, detectorPath + "/OyT_actual") * 1e-3;
  pos.shiftDown = m_loadHelper.getDoubleFromNexusPath(firstEntry, detectorPath + "/OyB_actual") * 1e-3;
236
  pos >> g_log.debug();
237
  return pos;
238
239
}

240
/**
241
 * Loads data for D11, D16 and D22
242
243
 * @param firstEntry : already opened first entry in nexus
 * @param instrumentPath : the path inside nexus where the instrument name is written
244
 */
Samuel Jones's avatar
Samuel Jones committed
245
void LoadILLSANS::initWorkSpace(NeXus::NXEntry &firstEntry, const std::string &instrumentPath) {
246
  g_log.debug("Fetching data...");
247
248
249
250
251
252
253
  std::string path;
  if (firstEntry.containsGroup("data")) {
    path = "data";
  } else {
    path = "data_scan/detector_data/data";
  }
  NXData dataGroup = firstEntry.openNXData(path);
254
255
  NXInt data = dataGroup.openIntData();
  data.load();
256
  size_t numberOfHistograms;
257

Samuel Jones's avatar
Samuel Jones committed
258
  m_isD16Omega = (data.dim0() == 1 && data.dim2() > 1 && m_instrumentName == "D16");
259
260

  if (m_isD16Omega) {
Samuel Jones's avatar
Samuel Jones committed
261
    numberOfHistograms = static_cast<size_t>(data.dim1() * data.dim2()) + N_MONITORS;
262
  } else {
Samuel Jones's avatar
Samuel Jones committed
263
    numberOfHistograms = static_cast<size_t>(data.dim0() * data.dim1()) + N_MONITORS;
264
  }
265
  createEmptyWorkspace(numberOfHistograms, 1);
266
  loadMetaData(firstEntry, instrumentPath);
267

268
  size_t nextIndex;
Mathieu Tillet's avatar
Mathieu Tillet committed
269
270
  nextIndex = loadDataFromTubes(data, m_defaultBinning, 0);
  nextIndex = loadDataFromMonitors(firstEntry, nextIndex);
271
  if (data.dim1() == 128) {
272
273
274
275
    m_resMode = "low";
  }
}

Mathieu Tillet's avatar
Mathieu Tillet committed
276
277
/**
 * @brief LoadILLSANS::initWorkSpaceD11B Load D11B data
278
279
 * @param firstEntry : already opened first entry in nexus
 * @param instrumentPath : the path inside nexus where the instrument name is written
Mathieu Tillet's avatar
Mathieu Tillet committed
280
 */
Samuel Jones's avatar
Samuel Jones committed
281
void LoadILLSANS::initWorkSpaceD11B(NeXus::NXEntry &firstEntry, const std::string &instrumentPath) {
Mathieu Tillet's avatar
Mathieu Tillet committed
282
283
  g_log.debug("Fetching data...");

284
  NXData data1 = firstEntry.openNXData("D11/Detector 1/data");
Mathieu Tillet's avatar
Mathieu Tillet committed
285
286
  NXInt dataCenter = data1.openIntData();
  dataCenter.load();
287
  NXData data2 = firstEntry.openNXData("D11/Detector 2/data");
288
  NXInt dataLeft = data2.openIntData();
Mathieu Tillet's avatar
Mathieu Tillet committed
289
  dataLeft.load();
290
  NXData data3 = firstEntry.openNXData("D11/Detector 3/data");
291
292
  NXInt dataRight = data3.openIntData();
  dataRight.load();
Mathieu Tillet's avatar
Mathieu Tillet committed
293
294

  size_t numberOfHistograms =
Samuel Jones's avatar
Samuel Jones committed
295
      static_cast<size_t>(dataCenter.dim0() * dataCenter.dim1() + dataRight.dim0() * dataRight.dim1() +
Mathieu Tillet's avatar
Mathieu Tillet committed
296
297
298
                          dataLeft.dim0() * dataLeft.dim1()) +
      N_MONITORS;

299
300
301
302
303
  if (dataCenter.dim2() != 1) {
    createEmptyWorkspace(numberOfHistograms, dataCenter.dim2(), MultichannelType::KINETIC);
  } else {
    createEmptyWorkspace(numberOfHistograms, 1);
  }
Mathieu Tillet's avatar
Mathieu Tillet committed
304
305
  loadMetaData(firstEntry, instrumentPath);

306
307
308
309
310
311
312
313
314
315
316
  // we need to adjust the default binning after loadmetadata
  if (dataCenter.dim2() != 1) {
    std::vector<double> frames(dataCenter.dim2(), 0);
    for (int i = 0; i < dataCenter.dim2(); ++i) {
      frames[i] = i;
    }
    m_defaultBinning.resize(dataCenter.dim2());
    std::copy(frames.cbegin(), frames.cend(), m_defaultBinning.begin());
  }

  MultichannelType type = (dataCenter.dim2() != 1) ? MultichannelType::KINETIC : MultichannelType::TOF;
Mathieu Tillet's avatar
Mathieu Tillet committed
317
  size_t nextIndex;
318
319
320
321
322
323
324
325
326
327
328
329
330
331
  nextIndex = loadDataFromTubes(dataCenter, m_defaultBinning, 0, type);
  nextIndex = loadDataFromTubes(dataLeft, m_defaultBinning, nextIndex, type);
  nextIndex = loadDataFromTubes(dataRight, m_defaultBinning, nextIndex, type);
  nextIndex = loadDataFromMonitors(firstEntry, nextIndex, type);
  if (dataCenter.dim2() != 1 && nextIndex < numberOfHistograms) {
    // there are a few runs with no 2nd monitor in kinetic, so we load the first monitor once again to preserve the
    // dimensions and x-axis
    nextIndex = loadDataFromMonitors(firstEntry, nextIndex, type);
  }
  // hijack the second monitor spectrum to store per-frame durations to enable time normalisation
  if (dataCenter.dim2() != 1) {
    NXFloat durations = firstEntry.openNXFloat("slices");
    durations.load();
    const HistogramData::Counts histoCounts(durations(), durations() + dataCenter.dim2());
332
    m_localWorkspace->setCounts(nextIndex - 1, histoCounts);
333
334
335
    m_localWorkspace->setCountVariances(nextIndex - 1,
                                        HistogramData::CountVariances(std::vector<double>(dataCenter.dim2(), 0)));
  }
Mathieu Tillet's avatar
Mathieu Tillet committed
336
337
338
339
}

/**
 * @brief LoadILLSANS::initWorkSpaceD22B Load D22B data
340
341
 * @param firstEntry : already opened first entry in nexus
 * @param instrumentPath : the path inside nexus where the instrument name is written
Mathieu Tillet's avatar
Mathieu Tillet committed
342
 */
Samuel Jones's avatar
Samuel Jones committed
343
void LoadILLSANS::initWorkSpaceD22B(NeXus::NXEntry &firstEntry, const std::string &instrumentPath) {
Mathieu Tillet's avatar
Mathieu Tillet committed
344
345
346
  g_log.debug("Fetching data...");

  NXData data2 = firstEntry.openNXData("data2");
347
348
  NXInt data2_data = data2.openIntData();
  data2_data.load();
349
  NXData data1 = firstEntry.openNXData("data1");
350
351
  NXInt data1_data = data1.openIntData();
  data1_data.load();
Mathieu Tillet's avatar
Mathieu Tillet committed
352
353

  size_t numberOfHistograms =
354
      static_cast<size_t>(data2_data.dim0() * data2_data.dim1() + data1_data.dim0() * data1_data.dim1()) + N_MONITORS;
Mathieu Tillet's avatar
Mathieu Tillet committed
355
356
357

  createEmptyWorkspace(numberOfHistograms, 1);
  loadMetaData(firstEntry, instrumentPath);
358
  runLoadInstrument();
Mathieu Tillet's avatar
Mathieu Tillet committed
359

360
  const std::string backIndex = m_localWorkspace->getInstrument()->getStringParameter("back_detector_index")[0];
Mathieu Tillet's avatar
Mathieu Tillet committed
361
  size_t nextIndex;
362
363
364
365
366
367
368
  if (backIndex == "2") {
    nextIndex = loadDataFromTubes(data2_data, m_defaultBinning, 0);
    nextIndex = loadDataFromTubes(data1_data, m_defaultBinning, nextIndex);
  } else {
    nextIndex = loadDataFromTubes(data1_data, m_defaultBinning, 0);
    nextIndex = loadDataFromTubes(data2_data, m_defaultBinning, nextIndex);
  }
Mathieu Tillet's avatar
Mathieu Tillet committed
369
  nextIndex = loadDataFromMonitors(firstEntry, nextIndex);
Mathieu Tillet's avatar
Mathieu Tillet committed
370
371
}

372
373
/**
 * Loads data for D33
374
375
 * @param firstEntry : already opened first entry in nexus
 * @param instrumentPath : the path inside nexus where the instrument name is written
376
 */
Samuel Jones's avatar
Samuel Jones committed
377
void LoadILLSANS::initWorkSpaceD33(NeXus::NXEntry &firstEntry, const std::string &instrumentPath) {
378
379
380
381
382
383
384

  g_log.debug("Fetching data...");

  NXData dataGroup1 = firstEntry.openNXData("data1");
  NXInt dataRear = dataGroup1.openIntData();
  dataRear.load();
  NXData dataGroup2 = firstEntry.openNXData("data2");
385
  NXInt dataRight = dataGroup2.openIntData();
386
  dataRight.load();
387
388
389
  NXData dataGroup3 = firstEntry.openNXData("data3");
  NXInt dataLeft = dataGroup3.openIntData();
  dataLeft.load();
390
391
392
393
394
395
396
397
398
  NXData dataGroup4 = firstEntry.openNXData("data4");
  NXInt dataDown = dataGroup4.openIntData();
  dataDown.load();
  NXData dataGroup5 = firstEntry.openNXData("data5");
  NXInt dataUp = dataGroup5.openIntData();
  dataUp.load();
  g_log.debug("Checking channel numbers...");

  // check number of channels
Samuel Jones's avatar
Samuel Jones committed
399
  if (dataRear.dim2() != dataRight.dim2() && dataRight.dim2() != dataLeft.dim2() &&
400
      dataLeft.dim2() != dataDown.dim2() && dataDown.dim2() != dataUp.dim2()) {
Samuel Jones's avatar
Samuel Jones committed
401
    throw std::runtime_error("The time bins have not the same dimension for all the 5 detectors!");
402
  }
403
  const auto numberOfHistograms = static_cast<size_t>(
Samuel Jones's avatar
Samuel Jones committed
404
405
      dataRear.dim0() * dataRear.dim1() + dataRight.dim0() * dataRight.dim1() + dataLeft.dim0() * dataLeft.dim1() +
      dataDown.dim0() * dataDown.dim1() + dataUp.dim0() * dataUp.dim1());
406
407

  g_log.debug("Creating empty workspace...");
Samuel Jones's avatar
Samuel Jones committed
408
  createEmptyWorkspace(numberOfHistograms + N_MONITORS, static_cast<size_t>(dataRear.dim2()));
409
410
411

  loadMetaData(firstEntry, instrumentPath);

Samuel Jones's avatar
Samuel Jones committed
412
  std::vector<double> binningRear, binningRight, binningLeft, binningDown, binningUp;
413
414
415
416
417
418
419
420
421

  if (firstEntry.getFloat("mode") == 0.0) { // Not TOF
    g_log.debug("Getting default wavelength bins...");
    binningRear = m_defaultBinning;
    binningRight = m_defaultBinning;
    binningLeft = m_defaultBinning;
    binningDown = m_defaultBinning;
    binningUp = m_defaultBinning;

Gagik Vardanyan's avatar
Gagik Vardanyan committed
422
423
  } else { // TOF
    m_isTOF = true;
Samuel Jones's avatar
Samuel Jones committed
424
    NXInt masterPair = firstEntry.openNXInt(m_instrumentName + "/tof/master_pair");
Gagik Vardanyan's avatar
Gagik Vardanyan committed
425
426
427
428
429
430
    masterPair.load();

    const std::string first = std::to_string(masterPair[0]);
    const std::string second = std::to_string(masterPair[1]);
    g_log.debug("Master choppers are " + first + " and " + second);

Samuel Jones's avatar
Samuel Jones committed
431
    NXFloat firstChopper = firstEntry.openNXFloat(m_instrumentName + "/chopper" + first + "/sample_distance");
Gagik Vardanyan's avatar
Gagik Vardanyan committed
432
    firstChopper.load();
Samuel Jones's avatar
Samuel Jones committed
433
    NXFloat secondChopper = firstEntry.openNXFloat(m_instrumentName + "/chopper" + second + "/sample_distance");
Gagik Vardanyan's avatar
Gagik Vardanyan committed
434
435
    secondChopper.load();
    m_sourcePos = (firstChopper[0] + secondChopper[0]) / 2.;
Samuel Jones's avatar
Samuel Jones committed
436
    g_log.debug("Source distance computed, moving moderator to Z=-" + std::to_string(m_sourcePos));
437
    g_log.debug("Getting wavelength bins from the nexus file...");
438
439
440
    bool vtof = true;
    // try VTOF mode
    try {
Samuel Jones's avatar
Samuel Jones committed
441
442
      NXInt channelWidthSum = firstEntry.openNXInt(m_instrumentName + "/tof/chwidth_sum");
      NXFloat channelWidthTimes = firstEntry.openNXFloat(m_instrumentName + "/tof/chwidth_times");
443
444
      channelWidthSum.load();
      channelWidthTimes.load();
445
      std::string distancePrefix(instrumentPath + "/tof/tof_distance_detector");
Samuel Jones's avatar
Samuel Jones committed
446
447
448
449
450
      binningRear = getVariableTimeBinning(firstEntry, distancePrefix + "1", channelWidthSum, channelWidthTimes);
      binningRight = getVariableTimeBinning(firstEntry, distancePrefix + "2", channelWidthSum, channelWidthTimes);
      binningLeft = getVariableTimeBinning(firstEntry, distancePrefix + "3", channelWidthSum, channelWidthTimes);
      binningDown = getVariableTimeBinning(firstEntry, distancePrefix + "4", channelWidthSum, channelWidthTimes);
      binningUp = getVariableTimeBinning(firstEntry, distancePrefix + "5", channelWidthSum, channelWidthTimes);
451
    } catch (const std::runtime_error &) {
452
453
454
455
456
      vtof = false;
    }
    if (!vtof) {
      try {
        // LTOF mode
Samuel Jones's avatar
Samuel Jones committed
457
458
459
460
461
462
        std::string binPathPrefix(instrumentPath + "/tof/tof_wavelength_detector");
        binningRear = m_loadHelper.getTimeBinningFromNexusPath(firstEntry, binPathPrefix + "1");
        binningRight = m_loadHelper.getTimeBinningFromNexusPath(firstEntry, binPathPrefix + "2");
        binningLeft = m_loadHelper.getTimeBinningFromNexusPath(firstEntry, binPathPrefix + "3");
        binningDown = m_loadHelper.getTimeBinningFromNexusPath(firstEntry, binPathPrefix + "4");
        binningUp = m_loadHelper.getTimeBinningFromNexusPath(firstEntry, binPathPrefix + "5");
463
      } catch (std::runtime_error &e) {
Samuel Jones's avatar
Samuel Jones committed
464
        throw std::runtime_error("Unable to load the wavelength axes for TOF data " + std::string(e.what()));
465
466
      }
    }
467
  }
468

469
  g_log.debug("Loading the data into the workspace...");
470

Mathieu Tillet's avatar
Mathieu Tillet committed
471
472
473
474
475
476
  size_t nextIndex = loadDataFromTubes(dataRear, binningRear, 0);
  nextIndex = loadDataFromTubes(dataRight, binningRight, nextIndex);
  nextIndex = loadDataFromTubes(dataLeft, binningLeft, nextIndex);
  nextIndex = loadDataFromTubes(dataDown, binningDown, nextIndex);
  nextIndex = loadDataFromTubes(dataUp, binningUp, nextIndex);
  nextIndex = loadDataFromMonitors(firstEntry, nextIndex);
477
478
}

479
480
481
482
483
484
485
486
/**
 * @brief Loads data from all the monitors
 * @param firstEntry : already opened first entry in nexus
 * @param firstIndex : the workspace index to load the first monitor to
 * @param type : used to discrimante between TOF and Kinetic
 * @return the next ws index after all the monitors
 */
size_t LoadILLSANS::loadDataFromMonitors(NeXus::NXEntry &firstEntry, size_t firstIndex, const MultichannelType type) {
487

488
  // let's find the monitors; should be monitor1 and monitor2
Samuel Jones's avatar
Samuel Jones committed
489
490
  for (std::vector<NXClassInfo>::const_iterator it = firstEntry.groups().begin(); it != firstEntry.groups().end();
       ++it) {
491
492
493
494
    if (it->nxclass == "NXmonitor") {
      NXData dataGroup = firstEntry.openNXData(it->nxname);
      NXInt data = dataGroup.openIntData();
      data.load();
Samuel Jones's avatar
Samuel Jones committed
495
496
      g_log.debug() << "Monitor: " << it->nxname << " dims = " << data.dim0() << "x" << data.dim1() << "x"
                    << data.dim2() << '\n';
497
      const HistogramData::Counts histoCounts(data(), data() + data.dim2());
498
      m_localWorkspace->setCounts(firstIndex, histoCounts);
499
      const HistogramData::CountVariances histoVariances(data(), data() + data.dim2());
500
      m_localWorkspace->setCountVariances(firstIndex, histoVariances);
501
502
503
504
505

      if (m_isTOF) {
        HistogramData::BinEdges histoBinEdges(data.dim2() + 1, HistogramData::LinearGenerator(0.0, 1));
        m_localWorkspace->setBinEdges(firstIndex, std::move(histoBinEdges));
      } else {
506
507
508
        if (type != MultichannelType::KINETIC) {
          HistogramData::BinEdges histoBinEdges = HistogramData::BinEdges(m_defaultBinning);
          m_localWorkspace->setBinEdges(firstIndex, std::move(histoBinEdges));
509
        } else {
510
511
          HistogramData::Points histoPoints = HistogramData::Points(m_defaultBinning);
          m_localWorkspace->setPoints(firstIndex, std::move(histoPoints));
512
        }
513
      }
514
      // Add average monitor counts to a property:
Samuel Jones's avatar
Samuel Jones committed
515
      double averageMonitorCounts = std::accumulate(data(), data() + data.dim2(), 0) / static_cast<double>(data.dim2());
516
517
518
519
520
521
522
523
524
      // make sure the monitor has values!
      if (averageMonitorCounts > 0) {
        API::Run &runDetails = m_localWorkspace->mutableRun();
        runDetails.addProperty("monitor", averageMonitorCounts, true);
      }
      firstIndex++;
    }
  }
  return firstIndex;
525
526
}

527
528
529
530
531
532
533
534
535
536
/**
 * @brief Loads data from tubes
 * @param data : a reference to already loaded nexus data block
 * @param timeBinning : the x-axis binning
 * @param firstIndex : the workspace index to start loading to
 * @param type : used to discrimante between TOF and Kinetic
 * @return the next ws index after all the tubes in the given detector bank
 */
size_t LoadILLSANS::loadDataFromTubes(NeXus::NXInt &data, const std::vector<double> &timeBinning, size_t firstIndex = 0,
                                      const MultichannelType type) {
537
  int numberOfTubes;
538
539
  int numberOfChannels;
  const int numberOfPixelsPerTube = data.dim1();
540

541
  if (m_isD16Omega) {
542
543
    // D16 with omega scan case
    numberOfTubes = data.dim2();
544
    numberOfChannels = data.dim0();
545
546
  } else {
    numberOfTubes = data.dim0();
547
    numberOfChannels = data.dim2();
548
  }
Mathieu Tillet's avatar
Mathieu Tillet committed
549

550
  PARALLEL_FOR_IF(Kernel::threadSafe(*m_localWorkspace))
551
552
  for (int i = 0; i < numberOfTubes; ++i) {
    for (int j = 0; j < numberOfPixelsPerTube; ++j) {
553
      int *data_p;
554
555
      if (m_isD16Omega) {
        data_p = &data(0, i, j);
556
      } else {
557
        data_p = &data(i, j, 0);
558
      }
559
      const size_t index = firstIndex + i * numberOfPixelsPerTube + j;
560
561
      const HistogramData::Counts histoCounts(data_p, data_p + numberOfChannels);
      const HistogramData::CountVariances histoVariances(data_p, data_p + numberOfChannels);
562
563
      m_localWorkspace->setCounts(index, histoCounts);
      m_localWorkspace->setCountVariances(index, histoVariances);
564

565
566
      if (type == MultichannelType::KINETIC) {
        const HistogramData::Points histoPoints(timeBinning);
567
        m_localWorkspace->setPoints(index, histoPoints);
568
      } else {
569
570
        const HistogramData::BinEdges binEdges(timeBinning);
        m_localWorkspace->setBinEdges(index, binEdges);
571
      }
572
573
    }
  }
574

575
  return firstIndex + numberOfTubes * numberOfPixelsPerTube;
576
577
}

578
/**
579
 * Create a workspace without any data in it
580
581
 * @param numberOfHistograms : number of spectra
 * @param numberOfChannels : number of TOF channels
582
 * @param type : type of the multichannel workspace (TOF (default) or Kinetic)
583
 */
584
585
586
void LoadILLSANS::createEmptyWorkspace(const size_t numberOfHistograms, const size_t numberOfChannels,
                                       const MultichannelType type) {
  const size_t numberOfElementsInX = numberOfChannels + ((type == MultichannelType::TOF && !m_isD16Omega) ? 1 : 0);
Samuel Jones's avatar
Samuel Jones committed
587
  m_localWorkspace =
588
589
590
591
      WorkspaceFactory::Instance().create("Workspace2D", numberOfHistograms, numberOfElementsInX, numberOfChannels);
  if (type == MultichannelType::TOF) {
    m_localWorkspace->getAxis(0)->unit() = UnitFactory::Instance().create("Wavelength");
  }
592
  m_localWorkspace->setYUnitLabel("Counts");
593
594
}

595
/**
LamarMoore's avatar
LamarMoore committed
596
597
598
599
600
 * Makes up the full path of the relevant IDF dependent on resolution mode
 * @param instName : the name of the instrument (including the resolution mode
 * suffix)
 * @return : the full path to the corresponding IDF
 */
Samuel Jones's avatar
Samuel Jones committed
601
std::string LoadILLSANS::getInstrumentFilePath(const std::string &instName) const {
602
603
604
605
606
607
608
609
610
611

  Poco::Path directory(ConfigService::Instance().getInstrumentDirectory());
  Poco::Path file(instName + "_Definition.xml");
  Poco::Path fullPath(directory, file);
  return fullPath.toString();
}

/**
 * Loads the instrument from the IDF
 */
612
613
void LoadILLSANS::runLoadInstrument() {

614
  auto loadInst = createChildAlgorithm("LoadInstrument");
615
  if (m_resMode == "nominal") {
Samuel Jones's avatar
Samuel Jones committed
616
    loadInst->setPropertyValue("Filename", getInstrumentFilePath(m_instrumentName));
617
  } else if (m_resMode == "low") {
618
    // low resolution mode we have only defined for the old D11 and D22
Samuel Jones's avatar
Samuel Jones committed
619
    loadInst->setPropertyValue("Filename", getInstrumentFilePath(m_instrumentName + "lr"));
620
  }
621
  loadInst->setProperty<MatrixWorkspace_sptr>("Workspace", m_localWorkspace);
Samuel Jones's avatar
Samuel Jones committed
622
  loadInst->setProperty("RewriteSpectraMap", Mantid::Kernel::OptionalBool(true));
623
  loadInst->execute();
624
625
}

626
627
628
629
/**
 * Move the detector banks for D33
 * @param detPos : structure holding the positions
 */
630
void LoadILLSANS::moveDetectorsD33(const DetectorPosition &detPos) {
631
632
633
634
635
636
637
638
639
640
641
642
  // Move in Z
  moveDetectorDistance(detPos.distanceSampleRear, "back_detector");
  moveDetectorDistance(detPos.distanceSampleBottomTop, "front_detector_top");
  moveDetectorDistance(detPos.distanceSampleBottomTop, "front_detector_bottom");
  moveDetectorDistance(detPos.distanceSampleRightLeft, "front_detector_right");
  moveDetectorDistance(detPos.distanceSampleRightLeft, "front_detector_left");
  // Move in X
  moveDetectorHorizontal(detPos.shiftLeft, "front_detector_left");
  moveDetectorHorizontal(-detPos.shiftRight, "front_detector_right");
  // Move in Y
  moveDetectorVertical(detPos.shiftUp, "front_detector_top");
  moveDetectorVertical(-detPos.shiftDown, "front_detector_bottom");
643
644
  // Set the sample log
  API::Run &runDetails = m_localWorkspace->mutableRun();
645
  runDetails.addProperty<double>("L2", detPos.distanceSampleRear, true);
646
647
648
649
}

/**
 * Move detectors in Z axis (X,Y are kept constant)
650
651
 * @param distance : the distance to move along Z axis [meters]
 * @param componentName : name of the component to move
652
 */
Samuel Jones's avatar
Samuel Jones committed
653
void LoadILLSANS::moveDetectorDistance(double distance, const std::string &componentName) {
654

655
  auto mover = createChildAlgorithm("MoveInstrumentComponent");
656
  V3D pos = getComponentPosition(componentName);
657
658
  mover->setProperty<MatrixWorkspace_sptr>("Workspace", m_localWorkspace);
  mover->setProperty("ComponentName", componentName);
659
660
661

  mover->setProperty("X", pos.X());
  mover->setProperty("Y", pos.Y());
662
  mover->setProperty("Z", distance);
663
  mover->setProperty("RelativePosition", false);
664
  mover->executeAsChildAlg();
665

Samuel Jones's avatar
Samuel Jones committed
666
  g_log.debug() << "Moving component '" << componentName << "' to Z = " << distance << '\n';
667
668
}

669
/**
670
 * Rotates instrument detector around y-axis in place
671
 * @param angle : the angle to rotate [degree]
672
 * @param componentName : "detector"
673
 */
Samuel Jones's avatar
Samuel Jones committed
674
void LoadILLSANS::rotateInstrument(double angle, const std::string &componentName) {
675
  auto rotater = createChildAlgorithm("RotateInstrumentComponent");
676
677
678
679
680
681
682
683
  rotater->setProperty<MatrixWorkspace_sptr>("Workspace", m_localWorkspace);
  rotater->setProperty("ComponentName", componentName);
  rotater->setProperty("X", 0.);
  rotater->setProperty("Y", 1.);
  rotater->setProperty("Z", 0.);
  rotater->setProperty("Angle", angle);
  rotater->setProperty("RelativeRotation", false);
  rotater->executeAsChildAlg();
Samuel Jones's avatar
Samuel Jones committed
684
  g_log.debug() << "Rotating component '" << componentName << "' to angle = " << angle << " degrees.\n";
685
686
}

687
688
689
690
691
692
/**
 * @brief LoadILLSANS::placeD16 : place the D16 detector.
 * @param angle : the angle between its center and the transmitted beam
 * @param distance : the distance between its center and the sample
 * @param componentName : "detector"
 */
Samuel Jones's avatar
Samuel Jones committed
693
void LoadILLSANS::placeD16(double angle, double distance, const std::string &componentName) {
694
  auto mover = createChildAlgorithm("MoveInstrumentComponent");
695
696
697
698
699
700
701
702
703
704
  mover->setProperty<MatrixWorkspace_sptr>("Workspace", m_localWorkspace);
  mover->setProperty("ComponentName", componentName);
  mover->setProperty("X", sin(angle * M_PI / 180) * distance);
  mover->setProperty("Y", 0.);
  mover->setProperty("Z", cos(angle * M_PI / 180) * distance);
  mover->setProperty("RelativePosition", false);
  mover->executeAsChildAlg();

  // rotate the detector so it faces the sample.
  rotateInstrument(angle, componentName);
705
706
  API::Run &runDetails = m_localWorkspace->mutableRun();
  runDetails.addProperty<double>("L2", distance, true);
707

Samuel Jones's avatar
Samuel Jones committed
708
  g_log.debug() << "Moving component '" << componentName << "' to angle = " << angle
709
710
711
                << " degrees and distance = " << distance << "metres.\n";
}

712
713
/**
 * Move detectors in X
714
715
 * @param shift : the distance to move [metres]
 * @param componentName : the name of the component
716
 */
Samuel Jones's avatar
Samuel Jones committed
717
void LoadILLSANS::moveDetectorHorizontal(double shift, const std::string &componentName) {
718
  auto mover = createChildAlgorithm("MoveInstrumentComponent");
719
  V3D pos = getComponentPosition(componentName);
720
721
  mover->setProperty<MatrixWorkspace_sptr>("Workspace", m_localWorkspace);
  mover->setProperty("ComponentName", componentName);
Gagik Vardanyan's avatar
Gagik Vardanyan committed
722
  mover->setProperty("X", shift);
723
724
725
726
  mover->setProperty("Y", pos.Y());
  mover->setProperty("Z", pos.Z());
  mover->setProperty("RelativePosition", false);
  mover->executeAsChildAlg();
Samuel Jones's avatar
Samuel Jones committed
727
  g_log.debug() << "Moving component '" << componentName << "' to X = " << shift << '\n';
728
729
}

730
731
732
733
734
/**
 * Move detectors in Y
 * @param shift : the distance to move [metres]
 * @param componentName : the name of the component
 */
Samuel Jones's avatar
Samuel Jones committed
735
void LoadILLSANS::moveDetectorVertical(double shift, const std::string &componentName) {
736
  auto mover = createChildAlgorithm("MoveInstrumentComponent");
737
  V3D pos = getComponentPosition(componentName);
738
739
740
741
742
743
744
  mover->setProperty<MatrixWorkspace_sptr>("Workspace", m_localWorkspace);
  mover->setProperty("ComponentName", componentName);
  mover->setProperty("X", pos.X());
  mover->setProperty("Y", shift);
  mover->setProperty("Z", pos.Z());
  mover->setProperty("RelativePosition", false);
  mover->executeAsChildAlg();
Samuel Jones's avatar
Samuel Jones committed
745
  g_log.debug() << "Moving component '" << componentName << "' to Y = " << shift << '\n';
746
747
748
}

/**
749
750
751
 * Get position of a component
 * @param componentName : the name of the component
 * @return : V3D of the component position
752
 */
753
V3D LoadILLSANS::getComponentPosition(const std::string &componentName) {
Samuel Jones's avatar
Samuel Jones committed
754
755
  Geometry::Instrument_const_sptr instrument = m_localWorkspace->getInstrument();
  Geometry::IComponent_const_sptr component = instrument->getComponentByName(componentName);
756
  return component->getPos();
757
758
}

759
760
761
/**
 * Loads some metadata present in the nexus file
 * @param entry : opened nexus entry
762
 * @param instrumentNamePath : the nexus entry of the instrument
763
 */
Samuel Jones's avatar
Samuel Jones committed
764
void LoadILLSANS::loadMetaData(const NeXus::NXEntry &entry, const std::string &instrumentNamePath) {
765
766
767
768

  g_log.debug("Loading metadata...");
  API::Run &runDetails = m_localWorkspace->mutableRun();

Samuel Jones's avatar
Samuel Jones committed
769
  if ((entry.getFloat("mode") == 0.0) || (m_instrumentName == "D16")) { // Not TOF
770
771
772
773
774
    runDetails.addProperty<std::string>("tof_mode", "Non TOF");
  } else {
    runDetails.addProperty<std::string>("tof_mode", "TOF");
  }

775
776
  double wavelength;
  if (getPointerToProperty("Wavelength")->isDefault()) {
777
778
779
780
781
    if (m_instrumentName == "D16") {
      wavelength = entry.getFloat(instrumentNamePath + "/Beam/wavelength");
    } else {
      wavelength = entry.getFloat(instrumentNamePath + "/selector/wavelength");
    }
Samuel Jones's avatar
Samuel Jones committed
782
    g_log.debug() << "Wavelength found in the nexus file: " << wavelength << '\n';
783
784
  } else {
    wavelength = getProperty("Wavelength");
785
  }
786

787
788
  // round the wavelength to avoid unnecessary rebinning during merge runs
  wavelength = std::round(wavelength * 100) / 100.;
789
790

  if (wavelength <= 0) {
791
    g_log.debug() << "Mode = " << entry.getFloat("mode") << '\n';
792
793
794
795
796
797
798
    g_log.information("The wavelength present in the NeXus file <= 0.");
    if (entry.getFloat("mode") == 0.0) { // Not TOF
      throw std::runtime_error("Working in Non TOF mode and the wavelength in "
                               "the file is <=0 !!! Check with the instrument "
                               "scientist!");
    }
  } else {
799
800
801
802
    double wavelengthRes = 10.;
    const std::string entryResolution = instrumentNamePath + "/selector/";
    try {
      wavelengthRes = entry.getFloat(entryResolution + "wavelength_res");
803
    } catch (const std::runtime_error &) {
804
805
      try {
        wavelengthRes = entry.getFloat(entryResolution + "wave_length_res");
806
      } catch (const std::runtime_error &) {
807
        if (m_instrumentName == "D16")
808
          wavelengthRes = 1;
Samuel Jones's avatar
Samuel Jones committed
809
        g_log.information() << "Could not find wavelength resolution, assuming " << wavelengthRes << "%.\n";
810
811
      }
    }
812
813
    // round also the wavelength res to avoid unnecessary rebinning during
    // merge runs
814
    wavelengthRes = std::round(wavelengthRes * 100) / 100.;
815
    runDetails.addProperty<double>("wavelength", wavelength);
Mathieu Tillet's avatar
Mathieu Tillet committed
816
    double ei = m_loadHelper.calculateEnergy(wavelength);
817
818
819
820
821
    runDetails.addProperty<double>("Ei", ei, true);
    // wavelength
    m_defaultBinning[0] = wavelength - wavelengthRes * wavelength * 0.01 / 2;
    m_defaultBinning[1] = wavelength + wavelengthRes * wavelength * 0.01 / 2;
  }
822
823
824
  // Add a log called timer with the value of duration
  const double duration = entry.getFloat("duration");
  runDetails.addProperty<double>("timer", duration);
825

826
  // the start time is needed in the workspace when loading the parameter file
827
  std::string startDate = entry.getString("start_time");
Samuel Jones's avatar
Samuel Jones committed
828
  runDetails.addProperty<std::string>("start_time", m_loadHelper.dateTimeInIsoFormat(startDate));
829
830
  // set the facility
  runDetails.addProperty<std::string>("Facility", std::string("ILL"));
831
832
}

833
834
835
836
/**
 * Sets full sample logs
 * @param filename : name of the file
 */
837
void LoadILLSANS::setFinalProperties(const std::string &filename) {
838
839
  API::Run &runDetails = m_localWorkspace->mutableRun();
  runDetails.addProperty("is_frame_skipping", 0);
840
841
  NXhandle nxHandle;
  NXstatus nxStat = NXopen(filename.c_str(), NXACC_READ, &nxHandle);
842

843
  if (nxStat != NX_ERROR) {
Mathieu Tillet's avatar
Mathieu Tillet committed
844
    m_loadHelper.addNexusFieldsToWsRun(nxHandle, runDetails);
845
    NXclose(&nxHandle);
846
  }
847
848
}

Gagik Vardanyan's avatar
Gagik Vardanyan committed
849
850
851
852
853
854
855
/**
 * Adjusts pixel by pixel the wavelength axis
 * Used only for D33 in TOF mode
 */
void LoadILLSANS::adjustTOF() {
  const auto &specInfo = m_localWorkspace->spectrumInfo();
  const double l1 = m_sourcePos;
856
  const size_t nHist = m_localWorkspace->getNumberHistograms();
Gagik Vardanyan's avatar
Gagik Vardanyan committed
857
  PARALLEL_FOR_IF(Kernel::threadSafe(*m_localWorkspace))
Samuel Jones's avatar
Samuel Jones committed
858
  for (int64_t index = 0; index < static_cast<int64_t>(nHist - N_MONITORS); ++index) {
<