Skip to content

Commit

Permalink
Refs #9564. Reversed direction of POLDI data.
Browse files Browse the repository at this point in the history
Unit tests had to be adjusted because of element ordering. While looking at the file,
made some performance improvements.
  • Loading branch information
Michael Wedel committed Jun 6, 2014
1 parent ff7108e commit 4eb1ba0
Show file tree
Hide file tree
Showing 4 changed files with 48 additions and 37 deletions.
Expand Up @@ -11,6 +11,7 @@
import mantid.simpleapi
from mantid import config
import os.path
import numpy as np

#--------- place to look for dictionary files

Expand Down Expand Up @@ -65,6 +66,16 @@ def PyExec(self):
ws = mantid.simpleapi.GroupDetectors(InputWorkspace=ws,
OutputWorkspace=wname,
MapFile=grp_file, Behaviour="Sum")

# Reverse direction of POLDI data so that low index corresponds to low 2theta.
histogramCount = ws.getNumberHistograms()
oldYData = []
for i in range(histogramCount):
oldYData.append([x for x in ws.readY(i)])

for i in range(histogramCount):
ws.setY(i, np.array(oldYData[histogramCount - 1 - i]))

elif inst == "TRICS":
ws = mantid.simpleapi.LoadFlexiNexus(fname,dicname,OutputWorkspace=wname)
ws = mantid.simpleapi.SINQTranspose3D(ws,OutputWorkspace=wname)
Expand Down
Expand Up @@ -49,36 +49,36 @@ class MANTID_SINQ_DLL PoldiAutoCorrelationCore
{
public:
PoldiAutoCorrelationCore(Kernel::Logger& g_log);
virtual ~PoldiAutoCorrelationCore() { }
~PoldiAutoCorrelationCore() { }

void setInstrument(boost::shared_ptr<PoldiAbstractDetector> detector, boost::shared_ptr<PoldiAbstractChopper> chopper);
void setWavelengthRange(double lambdaMin, double lambdaMax);

DataObjects::Workspace2D_sptr calculate(DataObjects::Workspace2D_sptr countData);

protected:
double getNormalizedTOFSum(std::vector<double> normalizedTofs);
std::vector<double> calculateDWeights(std::vector<double> tofsFor1Angstrom, double deltaT, double deltaD, size_t nd);
double getNormalizedTOFSum(const std::vector<double> &normalizedTofs) const;
std::vector<double> calculateDWeights(const std::vector<double> &tofsFor1Angstrom, double deltaT, double deltaD, size_t nd) const;

double getRawCorrelatedIntensity(double dValue, double weight);
virtual UncertainValue getCMessAndCSigma(double dValue, double slitTimeOffset, int index);
double reduceChopperSlitList(std::vector<UncertainValue> valuesWithSigma, double weight);
double getRawCorrelatedIntensity(double dValue, double weight) const;
UncertainValue getCMessAndCSigma(double dValue, double slitTimeOffset, int index) const;
double reduceChopperSlitList(const std::vector<UncertainValue> &valuesWithSigma, double weight) const;

std::vector<double> getDistances(std::vector<int> elements);
std::vector<double> getTofsFor1Angstrom(std::vector<int> elements);
std::vector<double> getDistances(const std::vector<int> &elements) const;
std::vector<double> getTofsFor1Angstrom(const std::vector<int> &elements) const;

virtual double getCounts(int x, int y);
virtual double getNormCounts(int x, int y);
double getCounts(int x, int y) const;
double getNormCounts(int x, int y) const;

virtual int getElementFromIndex(int index);
virtual double getTofFromIndex(int index);
double getSumOfCounts(int timeBinCount, std::vector<int> detectorElements);
int getElementFromIndex(int index) const;
double getTofFromIndex(int index) const;
double getSumOfCounts(int timeBinCount, const std::vector<int> &detectorElements) const;

virtual int cleanIndex(int index, int maximum);
int cleanIndex(int index, int maximum) const;

void setCountData(DataObjects::Workspace2D_sptr countData);

double correctedIntensity(double intensity, double weight);
double correctedIntensity(double intensity, double weight) const;


boost::shared_ptr<PoldiAbstractDetector> m_detector;
Expand Down
Expand Up @@ -77,7 +77,7 @@ DataObjects::Workspace2D_sptr PoldiAutoCorrelationCore::calculate(DataObjects::W
* - d-resolution deltaD, which results directly from deltaT
* - number of time bins for each copper cycle
*/
std::vector<double> timeData = m_countData->dataX(0);
std::vector<double> timeData = m_countData->readX(0);

m_logger.information() << " Setting time data..." << std::endl;
m_deltaT = timeData[1] - timeData[0];
Expand Down Expand Up @@ -182,7 +182,7 @@ DataObjects::Workspace2D_sptr PoldiAutoCorrelationCore::calculate(DataObjects::W
* @param normalizedTofs :: Vector with dimensionless numbers (TOF/(time bin size)).
* @return Sum over all elements of given vector.
*/
double PoldiAutoCorrelationCore::getNormalizedTOFSum(std::vector<double> normalizedTofs)
double PoldiAutoCorrelationCore::getNormalizedTOFSum(const std::vector<double> &normalizedTofs) const
{
/* Regarding numerical differences to fortran version
*
Expand All @@ -203,7 +203,7 @@ double PoldiAutoCorrelationCore::getNormalizedTOFSum(std::vector<double> normali
* @param nd :: Number of d-values.
* @return Vector with nd elements, each containing the weight for the given d-Value
*/
std::vector<double> PoldiAutoCorrelationCore::calculateDWeights(std::vector<double> tofsFor1Angstrom, double deltaT, double deltaD, size_t nd)
std::vector<double> PoldiAutoCorrelationCore::calculateDWeights(const std::vector<double> &tofsFor1Angstrom, double deltaT, double deltaD, size_t nd) const
{
/* Currently, all d-values get the same weight, so in fact this calculation would is not really
* necessary. But since there are not too many calculations involved and in order to stay close
Expand All @@ -227,7 +227,7 @@ std::vector<double> PoldiAutoCorrelationCore::calculateDWeights(std::vector<doub
* @param weight :: Dimensionless weight used for reduction of intermediate result.
* @return Correlated intensity, not corrected for correlation background.
*/
double PoldiAutoCorrelationCore::getRawCorrelatedIntensity(double dValue, double weight)
double PoldiAutoCorrelationCore::getRawCorrelatedIntensity(double dValue, double weight) const
{
/* Analysis according to the POLDI paper.
*
Expand Down Expand Up @@ -277,7 +277,7 @@ double PoldiAutoCorrelationCore::getRawCorrelatedIntensity(double dValue, double
* @param index :: Index of detector element the calculation is performed for.
* @return Pair of intensity and error for given input.
*/
UncertainValue PoldiAutoCorrelationCore::getCMessAndCSigma(double dValue, double slitTimeOffset, int index)
UncertainValue PoldiAutoCorrelationCore::getCMessAndCSigma(double dValue, double slitTimeOffset, int index) const
{
/* Element index and TOF for 1 Angstrom from current setup. */
int element = getElementFromIndex(index);
Expand Down Expand Up @@ -363,7 +363,7 @@ UncertainValue PoldiAutoCorrelationCore::getCMessAndCSigma(double dValue, double
* @param maximum :: Number of indices.
* @return Index on interval [0, max - 1]
*/
int PoldiAutoCorrelationCore::cleanIndex(int index, int maximum)
int PoldiAutoCorrelationCore::cleanIndex(int index, int maximum) const
{
int cleanIndex = index % maximum;
if(cleanIndex < 0) {
Expand All @@ -383,7 +383,7 @@ void PoldiAutoCorrelationCore::setCountData(DataObjects::Workspace2D_sptr countD
m_elementsMaxIndex = static_cast<int>(countData->getNumberHistograms()) - 1;
}

double PoldiAutoCorrelationCore::correctedIntensity(double intensity, double weight)
double PoldiAutoCorrelationCore::correctedIntensity(double intensity, double weight) const
{
return intensity - m_correlationBackground * weight / m_sumOfWeights;
}
Expand All @@ -394,7 +394,7 @@ double PoldiAutoCorrelationCore::correctedIntensity(double intensity, double wei
* @param weight :: Dimensionless weight.
* @return Correlated intensity, multiplied by weight.
*/
double PoldiAutoCorrelationCore::reduceChopperSlitList(std::vector<UncertainValue> valuesWithSigma, double weight)
double PoldiAutoCorrelationCore::reduceChopperSlitList(const std::vector<UncertainValue> &valuesWithSigma, double weight) const
{
try {
std::vector<double> signalToNoise(valuesWithSigma.size());
Expand All @@ -413,7 +413,7 @@ double PoldiAutoCorrelationCore::reduceChopperSlitList(std::vector<UncertainValu
* @param elements :: Vector of element indices.
* @return Vector of total neutron flight paths in mm.
*/
std::vector<double> PoldiAutoCorrelationCore::getDistances(std::vector<int> elements)
std::vector<double> PoldiAutoCorrelationCore::getDistances(const std::vector<int> &elements) const
{
double chopperDistance = m_chopper->distanceFromSample();
std::vector<double> distances;
Expand All @@ -433,7 +433,7 @@ std::vector<double> PoldiAutoCorrelationCore::getDistances(std::vector<int> elem
* @param elements :: Vector of element indices.
* @return Vector of specific TOFs (microseconds/Angstrom).
*/
std::vector<double> PoldiAutoCorrelationCore::getTofsFor1Angstrom(std::vector<int> elements)
std::vector<double> PoldiAutoCorrelationCore::getTofsFor1Angstrom(const std::vector<int> &elements) const
{
// Map element indices to 2Theta-Values
std::vector<double> twoThetas(elements.size());
Expand Down Expand Up @@ -462,9 +462,9 @@ std::vector<double> PoldiAutoCorrelationCore::getTofsFor1Angstrom(std::vector<in
* @param y :: Time bin index.
* @return Counts at position.
*/
double PoldiAutoCorrelationCore::getCounts(int x, int y)
double PoldiAutoCorrelationCore::getCounts(int x, int y) const
{
return static_cast<double>(m_countData->dataY(m_elementsMaxIndex - x)[y]);
return m_countData->readY(x)[y];
}

/** Returns normalized counts for correlation method at given position - these may come from a different source than the counts
Expand All @@ -473,9 +473,9 @@ double PoldiAutoCorrelationCore::getCounts(int x, int y)
* @param y :: Time bin index.
* @return Normalized counts at position.
*/
double PoldiAutoCorrelationCore::getNormCounts(int x, int y)
double PoldiAutoCorrelationCore::getNormCounts(int x, int y) const
{
return std::max(1.0, static_cast<double>(m_countData->dataY(m_elementsMaxIndex - x)[y]));
return std::max(1.0, m_countData->readY(x)[y]);
}

/** Returns detector element index for given index
Expand All @@ -485,7 +485,7 @@ double PoldiAutoCorrelationCore::getNormCounts(int x, int y)
* @param index :: Detector element index.
* @return Detector element index.
*/
int PoldiAutoCorrelationCore::getElementFromIndex(int index)
int PoldiAutoCorrelationCore::getElementFromIndex(int index) const
{
if(index < 0 || index >= static_cast<int>(m_detectorElements.size())) {
throw(std::range_error("Index out of bounds on accessing m_detectorElements."));
Expand All @@ -499,7 +499,7 @@ int PoldiAutoCorrelationCore::getElementFromIndex(int index)
* @param index :: Specific TOF index.
* @return TOF at index.
*/
double PoldiAutoCorrelationCore::getTofFromIndex(int index)
double PoldiAutoCorrelationCore::getTofFromIndex(int index) const
{
if(index < 0 || index >= static_cast<int>(m_tofsFor1Angstrom.size())) {
throw(std::range_error("Index out of bounds on accessing m_tofsFor1Angstrom."));
Expand All @@ -514,7 +514,7 @@ double PoldiAutoCorrelationCore::getTofFromIndex(int index)
* @param detectorElements :: Vector with all detector elements to be considered for the summation.
* @return Sum of all counts.
*/
double PoldiAutoCorrelationCore::getSumOfCounts(int timeBinCount, std::vector<int> detectorElements)
double PoldiAutoCorrelationCore::getSumOfCounts(int timeBinCount, const std::vector<int> &detectorElements) const
{
double sum = 0.0;

Expand Down
10 changes: 5 additions & 5 deletions Code/Mantid/Framework/SINQ/test/PoldiAutoCorrelationCoreTest.h
Expand Up @@ -214,10 +214,10 @@ class PoldiAutoCorrelationCoreTest : public CxxTest::TestSuite
TS_ASSERT_EQUALS(autoCorrelationCore.m_countData->getNumberHistograms(), 2);
TS_ASSERT_EQUALS(autoCorrelationCore.m_elementsMaxIndex, 1);

TS_ASSERT_EQUALS(autoCorrelationCore.getCounts(0, 0), 1.0);
TS_ASSERT_EQUALS(autoCorrelationCore.getCounts(0, 1), 1.0);
TS_ASSERT_EQUALS(autoCorrelationCore.getCounts(1, 0), 0.0);
TS_ASSERT_EQUALS(autoCorrelationCore.getCounts(1, 1), 0.0);
TS_ASSERT_EQUALS(autoCorrelationCore.getCounts(0, 0), 0.0);
TS_ASSERT_EQUALS(autoCorrelationCore.getCounts(0, 1), 0.0);
TS_ASSERT_EQUALS(autoCorrelationCore.getCounts(1, 0), 1.0);
TS_ASSERT_EQUALS(autoCorrelationCore.getCounts(1, 1), 1.0);
}

void testgetNormCounts()
Expand Down Expand Up @@ -263,7 +263,7 @@ class PoldiAutoCorrelationCoreTest : public CxxTest::TestSuite
double tofElements[] = {1.0, 2.0};
autoCorrelationCore.m_tofsFor1Angstrom = std::vector<double>(tofElements, tofElements + 2);

int elements[] = {1, 2};
int elements[] = {0, 1};
autoCorrelationCore.m_detectorElements = std::vector<int>(elements, elements + 2);

TS_ASSERT_DELTA(autoCorrelationCore.getCMessAndCSigma(1.2, 0.0, 0).value(), 0.0, 1e-6);
Expand Down

0 comments on commit 4eb1ba0

Please sign in to comment.