SANSCollimationLengthEstimator.cpp 7.12 KB
Newer Older
1
2
#include "MantidAlgorithms/SANSCollimationLengthEstimator.h"
#include "MantidAPI/MatrixWorkspace.h"
3
#include "MantidGeometry/Instrument.h"
4
5
#include "MantidKernel/Logger.h"
#include "MantidKernel/TimeSeriesProperty.h"
6
#include "MantidKernel/Property.h"
7
#include "MantidKernel/V3D.h"
8
#include "boost/lexical_cast.hpp"
9

10
11
namespace {
Mantid::Kernel::Logger g_log("SANSCollimationLengthEstimator");
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26

/**
 * Provide an string and check if it can be converted to a double
 * @param val: a value as a string
 * @returns true if it is convertible else false
 */
bool checkForDouble(std::string val) {
  auto isDouble = false;
  try {
    boost::lexical_cast<double>(val);
    isDouble = true;
  } catch (boost::bad_lexical_cast const &) {
  }
  return isDouble;
}
27
28
29
30
31
32
33
34
35
36
37
38
39
}

namespace Mantid {
namespace Algorithms {

using namespace Kernel;
using namespace API;

/**
 * Provide the collimation length which is associated with the instrument
 * @param workspace: the input workspace
 * @returns the collimation length
 */
40
41
42
43
double SANSCollimationLengthEstimator::provideCollimationLength(
    Mantid::API::MatrixWorkspace_sptr workspace) {
  // If the instrument does not have a correction specified then set the length
  // to 4
44
45
  const double defaultLColim = 4.0;
  auto collimationLengthID = "collimation-length-correction";
46

47
  if (!workspace->getInstrument()->hasParameter(collimationLengthID)) {
48
49
50
51
52
    g_log.error("Error in SANSCollimtionLengthEstimator: The instrument "
                "parameter file does not contain a collimation length "
                "correction,"
                "a default of 4 is provided. Please update the instrument "
                "parameter file.");
53
54
55
56
57
58
59
60
61
    return defaultLColim;
  }

  // Get the L1 length
  const V3D samplePos = workspace->getInstrument()->getSample()->getPos();
  const V3D sourcePos = workspace->getInstrument()->getSource()->getPos();
  const V3D SSD = samplePos - sourcePos;
  const double L1 = SSD.norm();

62
63
  auto collimationLengthCorrection =
      workspace->getInstrument()->getNumberParameter(collimationLengthID);
64

65
66
67
68
69
  if (workspace->getInstrument()->hasParameter(
          "special-default-collimation-length-method")) {
    auto specialCollimationMethod =
        workspace->getInstrument()->getStringParameter(
            "special-default-collimation-length-method");
70
    if (specialCollimationMethod[0] == "guide") {
71
72
73
74
75
76
77
78
79
      try {
        return getCollimationLengthWithGuides(workspace, L1,
                                              collimationLengthCorrection[0]);
      } catch (std::invalid_argument &ex) {
        g_log.notice() << ex.what();
        g_log.notice()
            << "SANSCollimationLengthEstimator: Not using any guides";
        return L1 - collimationLengthCorrection[0];
      }
80
    } else {
81
82
      throw std::invalid_argument("Error in SANSCollimationLengthEstimator: "
                                  "Unknown special collimation method.");
83
84
    }
  }
85
  return L1 - collimationLengthCorrection[0];
86
87
88
}

/**
89
90
91
92
 * This extraction strategy gets applied when guides are used to calculate the
 * collimation length. The instrument
 * parameter file contains the information about the number of guides to use.
 * The guide data itself is fetched
93
 * from the sample logs.
94
95
 * @param inOutWS: a matrix workspace
 * @param L1: the distance between sample and source
96
97
 * @param collimationLengthCorrection: The correction to get the collimation
 * length
98
 */
99
100
101
double SANSCollimationLengthEstimator::getCollimationLengthWithGuides(
    MatrixWorkspace_sptr inOutWS, const double L1,
    const double collimationLengthCorrection) const {
102
103
104
105
  auto lCollim = L1 - collimationLengthCorrection;

  // Make sure we have guide cutoffs
  if (!inOutWS->getInstrument()->hasParameter("guide-cutoff")) {
106
107
    throw std::invalid_argument("TOFSANSResolutionByPixel: Could not get a "
                                "GuideCutoff from the instrument");
108
109
110
111
  }

  // Make sure we have a defined number of guidess
  if (!inOutWS->getInstrument()->hasParameter("number-of-guides")) {
112
113
    throw std::invalid_argument(
        "TOFSANSResolutionByPixel: Could not get the number of guides.");
114
115
  }

116
117
118
119
  // Make sure we have a guide increment specified
  if (!inOutWS->getInstrument()->hasParameter(
          "guide-collimation-length-increment")) {
    throw std::invalid_argument(
120
        "TOFSANSResolutionByPixel: Could not find a guide increment.");
121
122
  }

123
124
125
126
  auto numberOfGuides = static_cast<unsigned int>(
      inOutWS->getInstrument()->getNumberParameter("number-of-guides")[0]);
  auto guideIncrement = inOutWS->getInstrument()->getNumberParameter(
      "guide-collimation-length-increment")[0];
127

128
129
  // Make sure that all guides are there. They are labelled as Guide1, Guide2,
  // Guide3, ...
130
131
  // The entry is a numeric TimeSeriesProperty or a numeric entry, if something
  // else then default
132
133
  std::vector<double> guideValues;
  for (unsigned int i = 1; i <= numberOfGuides; i++) {
134
    auto guideName = "Guide" + std::to_string(i);
135
    if (inOutWS->run().hasProperty(guideName)) {
136
137
      auto guideValue = getGuideValue(inOutWS->run().getProperty(guideName));
      guideValues.push_back(guideValue);
138
    } else {
139
140
141
      throw std::invalid_argument("TOFSANSResolutionByPixel: Mismatch between "
                                  "specified number of Guides and actual "
                                  "Guides.");
142
143
144
    }
  }

145
146
147
148
149
150
151
152
  auto guideCutoff =
      inOutWS->getInstrument()->getNumberParameter("guide-cutoff")[0];
  // Go through the guides and check in an alternate manner if the guide is
  // smaller
  // or larger than the cut off value. We start at the last guide and check that
  // it is
  // larger than the cutoff, the next one has to be smaller and so on. For
  // example in pseudocode
153
154
155
156
157
158
  // If Guide5 > 130: LCollim+=2.0 else break;
  // If Guide4 < 130: LCollim+=2.0 else break;
  // If Guide3 > 130: LCollim+=2.0 else break;
  // ...
  unsigned int largerSmallerCounter = 0;
  for (auto it = guideValues.rbegin(); it != guideValues.rend(); ++it) {
159
    bool guideIsLarger = largerSmallerCounter % 2 == 0;
160
    if (guideIsLarger && (*it > guideCutoff)) {
161
      lCollim += guideIncrement;
162
    } else if (!guideIsLarger && (*it < guideCutoff)) {
163
      lCollim += guideIncrement;
164
165
166
167
168
169
170
    } else {
      break;
    }
    largerSmallerCounter++;
  }
  return lCollim;
}
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195

/**
 * Extracts the value of the guide
 * @param prop: a property
 * @returns the guide value
 */
double SANSCollimationLengthEstimator::getGuideValue(
    Mantid::Kernel::Property *prop) const {
  if (auto timeSeriesProperty =
          dynamic_cast<TimeSeriesProperty<double> *>(prop)) {
    return timeSeriesProperty->firstValue();
  } else if (auto doubleProperty =
                 dynamic_cast<PropertyWithValue<double> *>(prop)) {
    auto val = doubleProperty->value();
    if (checkForDouble(val)) {
      g_log.warning("SANSCollimationLengthEstimator: The Guide was not "
                    "recoginized as a TimeSeriesProperty, but rather as a "
                    "Numeric.");
      return boost::lexical_cast<double, std::string>(val);
    }
  }
  throw std::invalid_argument("TOFSANSResolutionByPixel: Unknown type for "
                              "Guides. Currently only Numeric and TimeSeries "
                              "are supported.");
}
196
197
}
}