30using namespace Kernel;
32using namespace Geometry;
33using namespace DataObjects;
52 return "Bins TOF powder diffraction event data in 2D space.";
59 auto wsValidator = std::make_shared<CompositeValidator>();
67 "An input EventWorkspace must be a Histogram workspace, not Point data. "
68 "X-axis units must be wavelength.");
71 "An output workspace.");
73 const std::string docString =
"A comma separated list of first bin boundary, width, last bin boundary. "
75 "this can be followed by a comma and more widths and last boundary "
77 "Negative width values indicate logarithmic binning.";
78 auto rebinValidator = std::make_shared<RebinParamsValidator>(
true);
82 const std::vector<std::string> exts{
".txt",
".dat"};
84 "Optional: The ascii file containing the list of bin edges. "
85 "Either this or Axis1- and dPerpendicularBinning need to be specified.");
88 "Normalize the binned workspace by the bin area.");
101 const bool normalizeByBinArea = this->
getProperty(
"NormalizeByBinArea");
102 if (normalizeByBinArea)
113 std::map<std::string, std::string> result;
118 if (!useBinFile && !useBinning1 && !useBinning2) {
119 const std::string msg =
"You must specify either dSpaceBinning and "
120 "dPerpendicularBinning, or a BinEdgesFile.";
121 result[
"dSpaceBinning"] = msg;
122 result[
"dPerpendicularBinning"] = msg;
123 result[
"BinEdgesFile"] = msg;
124 }
else if (useBinFile && (useBinning1 || useBinning2)) {
125 const std::string msg =
"You must specify either dSpaceBinning and "
126 "dPerpendicularBinning, or a BinEdgesFile, but not both.";
127 result[
"BinEdgesFile"] = msg;
142 bool binsFromFile(
false);
143 size_t dPerpSize = 0;
146 const auto &spectrumInfo =
m_inputWS->spectrumInfo();
148 const std::string beFileName =
getProperty(
"BinEdgesFile");
149 if (!beFileName.empty())
153 BinEdges dBins(oldXEdges.size());
154 BinEdges dPerpBins(oldXEdges.size());
156 auto &dPerp = dPerpBins.mutableRawData();
157 std::vector<std::vector<double>> fileXbins;
163 dPerpSize = dPerp.size();
166 g_log.
debug() <<
"Maximal size of Xbins = " << dSize;
168 g_log.
debug() <<
"Outws has " << outputWS->getNumberHistograms() <<
" histograms and " << outputWS->blocksize()
169 <<
" bins." << std::endl;
172 for (
const auto &Xbins : fileXbins) {
173 g_log.
debug() <<
"Xbins size: " << Xbins.size() << std::endl;
174 BinEdges binEdges(Xbins);
175 outputWS->setBinEdges(idx, binEdges);
181 HistogramData::BinEdges binEdges(dBins);
183 dSize = binEdges.size();
185 for (
size_t idx = 0; idx < dPerpSize - 1; idx++)
186 outputWS->setBinEdges(idx, binEdges);
187 auto abscissa = std::make_unique<BinEdgeAxis>(dBins.mutableRawData());
188 outputWS->replaceAxis(0, std::move(abscissa));
193 auto verticalAxis = std::make_unique<BinEdgeAxis>(dPerp);
194 auto verticalAxisRaw = verticalAxis.get();
197 verticalAxis->title() =
"d_p";
198 outputWS->replaceAxis(1, std::move(verticalAxis));
202 std::vector<std::vector<double>> newYValues(dPerpSize - 1, std::vector<double>(dSize - 1, 0.0));
203 std::vector<std::vector<double>> newEValues(dPerpSize - 1, std::vector<double>(dSize - 1, 0.0));
206 g_log.
debug() <<
"newYSize = " << dPerpSize << std::endl;
207 g_log.
debug() <<
"newXSize = " << dSize << std::endl;
208 std::vector<double> dp_vec(verticalAxisRaw->getValues());
211 for (int64_t snum = 0; snum < numSpectra; ++snum) {
213 if (!spectrumInfo.isMasked(snum)) {
214 double theta = 0.5 * spectrumInfo.twoTheta(snum);
215 double sin_theta = sin(theta);
216 if (sin_theta == 0) {
217 throw std::runtime_error(
"Spectrum " +
std::to_string(snum) +
" has sin(theta)=0. Cannot calculate d-Spacing!");
219 if (cos(theta) <= 0) {
221 " has cos(theta) <= 0. Cannot calculate d-SpacingPerpendicular!");
223 double log_cos_theta = log(cos(theta));
232 for (
const auto &ev : events) {
233 auto d =
calcD(ev.tof(), sin_theta);
234 auto dp =
calcDPerp(ev.tof(), log_cos_theta);
235 const auto lowy = std::lower_bound(dp_vec.begin(), dp_vec.end(), dp);
236 if ((lowy == dp_vec.end()) || (lowy == dp_vec.begin()))
238 int64_t dp_index = std::distance(dp_vec.begin(), lowy) - 1;
239 auto xs = binsFromFile ? fileXbins[dp_index] : dBins.rawData();
240 const auto lowx = std::lower_bound(xs.begin(), xs.end(),
d);
241 if ((lowx == xs.end()) || lowx == xs.begin())
243 int64_t d_index = std::distance(xs.begin(), lowx) - 1;
246 newYValues[dp_index][d_index] += ev.weight();
247 newEValues[dp_index][d_index] += ev.errorSquared();
251 prog.
report(
"Binning event data...");
256 for (
const auto &yVec : newYValues) {
257 outputWS->setCounts(idx, yVec);
261 for (
auto &eVec : newEValues) {
262 std::transform(eVec.begin(), eVec.end(), eVec.begin(),
static_cast<double (*)(
double)
>(sqrt));
263 outputWS->setCountStandardDeviations(idx, eVec);
277 std::vector<std::vector<double>> &Xbins)
const {
278 const std::string beFileName =
getProperty(
"BinEdgesFile");
279 std::ifstream file(beFileName);
281 std::string::size_type
n;
282 std::string::size_type sz;
283 std::vector<double>
tmp;
285 while (getline(file, line)) {
286 n = line.find(
"dp =");
287 if (
n != std::string::npos) {
289 Xbins.emplace_back(
tmp);
292 double dp1 = std::stod(line.substr(4), &sz);
293 double dp2 = std::stod(line.substr(sz + 4));
295 Ybins.emplace_back(dp1);
296 Ybins.emplace_back(dp2);
298 Ybins.emplace_back(dp2);
301 }
else if (line.find(
"#") == std::string::npos) {
302 std::stringstream ss(line);
309 Xbins.emplace_back(
tmp);
311 g_log.
information() <<
"Number of Xbins sets: " << Xbins.size() << std::endl;
326 for (
const auto &v : Xbins) {
327 max_size = std::max(v.size(), max_size);
330 for (
auto &v : Xbins) {
331 if (v.size() < max_size)
332 v.resize(max_size, v.back());
339 const std::vector<double> &yValues = verticalAxis->
getValues();
340 auto nhist = outWS->getNumberHistograms();
341 g_log.
debug() <<
"Number of hists: " << nhist <<
" Length of YAxis: " << verticalAxis->
length() << std::endl;
343 for (
size_t idx = 0; idx < nhist; ++idx) {
344 double factor = 1.0 / (yValues[idx + 1] - yValues[idx]);
346 outWS->convertToFrequencies(idx);
347 auto &freqs = outWS->mutableY(idx);
348 using std::placeholders::_1;
349 std::transform(freqs.begin(), freqs.end(), freqs.begin(), std::bind(std::multiplies<double>(), factor, _1));
350 auto &errors = outWS->mutableE(idx);
351 std::transform(errors.begin(), errors.end(), errors.begin(), std::bind(std::multiplies<double>(), factor, _1));
355double calcD(
double wavelength,
double sintheta) {
return wavelength * 0.5 / sintheta; }
357double calcDPerp(
double wavelength,
double logcostheta) {
return sqrt(wavelength * wavelength - 2.0 * logcostheta); }
#define DECLARE_ALGORITHM(classname)
#define PARALLEL_START_INTERRUPT_REGION
Begins a block to skip processing is the algorithm has been interupted Note the end of the block if n...
#define PARALLEL_CRITICAL(name)
#define PARALLEL_END_INTERRUPT_REGION
Ends a block to skip processing is the algorithm has been interupted Note the start of the block if n...
#define PARALLEL_FOR_IF(condition)
Empty definitions - to enable set your complier to enable openMP.
#define PARALLEL_CHECK_INTERRUPT_REGION
Adds a check after a Parallel region to see if it was interupted.
void declareProperty(std::unique_ptr< Kernel::Property > p, const std::string &doc="") override
Add a property to the list of managed properties.
Kernel::Property * getPointerToProperty(const std::string &name) const override
Get a property by name.
TypedValue getProperty(const std::string &name) const override
Get the value of a property.
@ OptionalLoad
to specify a file to read but the file doesn't have to exist
A validator which checks that a workspace contains histogram data (the default) or point data as requ...
A validator which checks that a workspace has a valid instrument.
Class to represent a numeric axis of a workspace.
virtual const std::vector< double > & getValues() const
Return a const reference to the values.
std::size_t length() const override
Get the length of the axis.
Helper class for reporting progress from algorithms.
A validator which checks whether the input workspace has the Spectra number in the axis.
A property class for workspaces.
A validator which checks that the unit of the workspace referred to by a WorkspaceProperty is the exp...
void ReadBinsFromFile(std::vector< double > &Ybins, std::vector< std::vector< double > > &Xbins) const
Bin2DPowderDiffraction::ReadBinsFromFile.
const std::string summary() const override
Algorithm's summary for use in the GUI and help.
void normalizeToBinArea(const API::MatrixWorkspace_sptr &outWS)
void exec() override
Execute the algorithm.
int m_numberOfSpectra
The number of spectra in the workspace.
std::map< std::string, std::string > validateInputs() override
Cross-check properties with each other.
API::MatrixWorkspace_sptr createOutputWorkspace()
Setup the output workspace.
void init() override
Initialize the algorithm's properties.
int version() const override
Algorithm's version for identification.
size_t UnifyXBins(std::vector< std::vector< double > > &Xbins) const
Bin2DPowderDiffraction::UnifyXBins unifies size of the vectors in Xbins.
DataObjects::EventWorkspace_sptr m_inputWS
Pointer to the input event workspace.
const std::string category() const override
Algorithm's category for identification.
Mantid::API::EventType getEventType() const override
Return the type of Event vector contained within.
void switchTo(Mantid::API::EventType newType) override
Switch the EventList to use the given EventType (TOF, WEIGHTED, or WEIGHTED_NOTIME)
std::vector< WeightedEvent > & getWeightedEvents()
Return the list of WeightedEvent contained.
Support for a property that holds an array of values.
IPropertyManager * setProperty(const std::string &name, const T &value)
Templated method to set the value of a PropertyWithValue.
void debug(const std::string &msg)
Logs at debug level.
void information(const std::string &msg)
Logs at information level.
void report()
Increments the loop counter by 1, then sends the progress notification on behalf of its algorithm.
The concrete, templated class for properties.
virtual bool isDefault() const =0
Overriden function that returns if property has the same value that it was initialised with,...
static T & Instance()
Return a reference to the Singleton instance, creating it if it does not already exist Creation is do...
std::shared_ptr< MatrixWorkspace > MatrixWorkspace_sptr
shared pointer to the matrix workspace base class
double calcDPerp(double wavelength, double logcostheta)
double calcD(double wavelength, double sintheta)
int MANTID_KERNEL_DLL createAxisFromRebinParams(const std::vector< double > ¶ms, std::vector< double > &xnew, const bool resize_xnew=true, const bool full_bins_only=false, const double xMinHint=std::nan(""), const double xMaxHint=std::nan(""), const bool useReverseLogarithmic=false, const double power=-1)
Creates a new output X array given a 'standard' set of rebinning parameters.
std::enable_if< std::is_pointer< Arg >::value, bool >::type threadSafe(Arg workspace)
Thread-safety check Checks the workspace to ensure it is suitable for multithreaded access.
std::string to_string(const wide_integer< Bits, Signed > &n)
@ Input
An input workspace.
@ Output
An output workspace.