Commit 1f5c8337 authored by Owen Arnold's avatar Owen Arnold
Browse files

refs #11056. New method createPeakHKL.

Unit tested to ensure we are creating a self-consistent peak.
Expose to python and unit test there too.
parent 941ce153
......@@ -51,7 +51,7 @@ public:
virtual bool findDetector() = 0;
virtual void setQSampleFrame(Mantid::Kernel::V3D QSampleFrame,
double detectorDistance = 1.0) = 0;
boost::optional<double> detectorDistance = boost::optional<double>()) = 0;
virtual void setQLabFrame(Mantid::Kernel::V3D QLabFrame,
boost::optional<double> detectorDistance = boost::optional<double>()) = 0;
......
......@@ -7,6 +7,7 @@
#include "MantidAPI/ITableWorkspace.h"
#include "MantidAPI/ExperimentInfo.h"
#include "MantidKernel/SpecialCoordinateSystem.h"
#include <boost/optional.hpp>
namespace Mantid {
......@@ -98,12 +99,20 @@ public:
//---------------------------------------------------------------------------------------------
/** Create an instance of a Peak
* @param QLabFrame :: Q of the center of the peak, in reciprocal space
* @param detectorDistance :: distance between the sample and the detector.
* @param QLabFrame :: Q of the center of the peak in the lab frame, in reciprocal space
* @param detectorDistance :: Optional distance between the sample and the detector. Calculated if not provided.
* @return a pointer to a new Peak object.
*/
virtual IPeak *createPeak(Mantid::Kernel::V3D QLabFrame,
double detectorDistance = 1.0) const = 0;
boost::optional<double> detectorDistance) const = 0;
/**
* Create an instance of a peak using a V3D
* @param HKL V3D
* @return a pointer to a new Peak object.
*/
virtual IPeak *createPeakHKL(Mantid::Kernel::V3D HKL) const = 0;
//---------------------------------------------------------------------------------------------
/** Determine if the workspace has been integrated using a peaks integration
......
......@@ -28,10 +28,10 @@ public:
Peak();
Peak(Geometry::Instrument_const_sptr m_inst, Mantid::Kernel::V3D QLabFrame,
double detectorDistance = 1.0);
boost::optional<double> detectorDistance = boost::optional<double>());
Peak(Geometry::Instrument_const_sptr m_inst, Mantid::Kernel::V3D QSampleFrame,
Mantid::Kernel::Matrix<double> goniometer,
double detectorDistance = 1.0);
boost::optional<double> detectorDistance = boost::optional<double>());
Peak(Geometry::Instrument_const_sptr m_inst, int m_DetectorID,
double m_Wavelength);
Peak(Geometry::Instrument_const_sptr m_inst, int m_DetectorID,
......@@ -86,7 +86,7 @@ public:
Mantid::Kernel::V3D getDetectorPositionNoCheck() const;
void setQSampleFrame(Mantid::Kernel::V3D QSampleFrame,
double detectorDistance = 1.0);
boost::optional<double> detectorDistance = boost::optional<double>());
void setQLabFrame(Mantid::Kernel::V3D QLabFrame,
boost::optional<double> detectorDistance = boost::optional<double>());
......
......@@ -21,6 +21,7 @@
#include <string>
#include <utility>
#include <vector>
#include <boost/optional.hpp>
// IsamplePosition should be IsampleOrientation
namespace Mantid {
......@@ -108,9 +109,13 @@ public:
const Peak &getPeak(int peakNum) const;
API::IPeak *createPeak(Kernel::V3D QFrame,
double detectorDistance = 1.0) const;
boost::optional<double> detectorDistance = boost::optional<double>()) const;
std::vector<std::pair<std::string, std::string>>
peakInfo(Kernel::V3D qFrame, bool labCoords) const;
Peak *createPeakHKL(Kernel::V3D HKL) const;
int peakInfoNumber(Kernel::V3D qFrame, bool labCoords) const;
std::vector<Peak> &getPeaks();
......
......@@ -30,11 +30,11 @@ Peak::Peak()
*
* @param m_inst :: Shared pointer to the instrument for this peak detection
* @param QLabFrame :: Q of the center of the peak, in reciprocal space
* @param detectorDistance :: distance between the sample and the detector.
* @param detectorDistance :: Optional distance between the sample and the detector. Calculated if not explicitly provided.
* Used to give a valid TOF. Default 1.0 meters.
*/
Peak::Peak(Geometry::Instrument_const_sptr m_inst,
Mantid::Kernel::V3D QLabFrame, double detectorDistance)
Mantid::Kernel::V3D QLabFrame, boost::optional<double> detectorDistance)
: m_H(0), m_K(0), m_L(0), m_Intensity(0), m_SigmaIntensity(0),
m_BinCount(0), m_GoniometerMatrix(3, 3, true),
m_InverseGoniometerMatrix(3, 3, true), m_RunNumber(0), m_MonitorCount(0),
......@@ -52,12 +52,12 @@ Peak::Peak(Geometry::Instrument_const_sptr m_inst,
* @param QSampleFrame :: Q of the center of the peak, in reciprocal space, in
*the sample frame (goniometer rotation accounted for).
* @param goniometer :: a 3x3 rotation matrix
* @param detectorDistance :: distance between the sample and the detector.
* @param detectorDistance :: Optional distance between the sample and the detector. Calculated if not explicitly provided.
* Used to give a valid TOF. Default 1.0 meters.
*/
Peak::Peak(Geometry::Instrument_const_sptr m_inst,
Mantid::Kernel::V3D QSampleFrame,
Mantid::Kernel::Matrix<double> goniometer, double detectorDistance)
Mantid::Kernel::Matrix<double> goniometer, boost::optional<double> detectorDistance)
: m_H(0), m_K(0), m_L(0), m_Intensity(0), m_SigmaIntensity(0),
m_BinCount(0), m_GoniometerMatrix(goniometer),
m_InverseGoniometerMatrix(goniometer), m_RunNumber(0), m_MonitorCount(0),
......@@ -464,10 +464,10 @@ Mantid::Kernel::V3D Peak::getQSampleFrame() const {
* This is in inelastic convention: momentum transfer of the LATTICE!
* Also, q does NOT have a 2pi factor = it is equal to 1/wavelength.
* @param detectorDistance :: distance between the sample and the detector.
* Used to give a valid TOF. Default 1.0 meters.
* Used to give a valid TOF. You do NOT need to explicitly set this.
*/
void Peak::setQSampleFrame(Mantid::Kernel::V3D QSampleFrame,
double detectorDistance) {
boost::optional<double> detectorDistance) {
V3D Qlab = m_GoniometerMatrix * QSampleFrame;
this->setQLabFrame(Qlab, detectorDistance);
}
......
......@@ -199,11 +199,11 @@ const Peak &PeaksWorkspace::getPeak(const int peakNum) const {
//---------------------------------------------------------------------------------------------
/** Creates an instance of a Peak BUT DOES NOT ADD IT TO THE WORKSPACE
* @param QLabFrame :: Q of the center of the peak, in reciprocal space
* @param detectorDistance :: distance between the sample and the detector.
* @param detectorDistance :: optional distance between the sample and the detector. You do NOT need to explicitly provide this distance.
* @return a pointer to a new Peak object.
*/
API::IPeak *PeaksWorkspace::createPeak(Kernel::V3D QLabFrame,
double detectorDistance) const {
boost::optional<double> detectorDistance) const {
return new Peak(this->getInstrument(), QLabFrame, detectorDistance);
}
......@@ -299,6 +299,7 @@ PeaksWorkspace::peakInfo(Kernel::V3D qFrame, bool labCoords) const {
}
try {
API::IPeak *peak = createPeak(Qlab);
if (sample().hasOrientedLattice()) {
......@@ -387,6 +388,33 @@ PeaksWorkspace::peakInfo(Kernel::V3D qFrame, bool labCoords) const {
return Result;
}
/**
* Create a Peak from a HKL value provided by the client.
*
*
* @param HKL : reciprocal lattice vector coefficients
* @return Fully formed peak.
*/
Peak *PeaksWorkspace::createPeakHKL(V3D HKL) const
{
Geometry::OrientedLattice lattice = this->sample().getOrientedLattice();
Geometry::Goniometer goniometer = this->run().getGoniometer();
// Calculate qLab from q HKL. As per Busing and Levy 1967, q_lab_frame = 2pi * Goniometer * UB * HKL
V3D qLabFrame = goniometer.getR() * lattice.getUB() * HKL * 2 * M_PI;
// create a peak using the qLab frame
auto peak = new Peak(this->getInstrument(), qLabFrame); // This should calculate the detector positions too.
// We need to set HKL separately to keep things consistent.
peak->setHKL(HKL[0], HKL[1], HKL[2]);
// Set the goniometer
peak->setGoniometerMatrix(goniometer.getR());
return peak;
}
/**
* Returns selected information for a "peak" at QLabFrame.
*
......
......@@ -271,8 +271,9 @@ public:
void test_setQLabFrame_ThrowsIfQIsNull()
{
Peak p1(inst, 10000, 2.0);
TS_ASSERT_THROWS_ANYTHING(Peak p2(inst, V3D(0,0,0), 1.0));
TS_ASSERT_THROWS_ANYTHING(Peak p2(inst, V3D(1,2,0), 1.0));
const boost::optional<double> distance = 1.0;
TS_ASSERT_THROWS_ANYTHING(Peak p2(inst, V3D(0,0,0), distance));
TS_ASSERT_THROWS_ANYTHING(Peak p2(inst, V3D(1,2,0), distance));
}
......@@ -301,7 +302,7 @@ public:
V3D detPos1 = p1.getDetPos();
// Construct using just Q
Peak p2(inst, Qlab1, detPos1.norm());
Peak p2(inst, Qlab1, boost::optional<double>(detPos1.norm()));
comparePeaks(p1, p2);
TS_ASSERT_EQUALS( p2.getBankName(), "None");
TS_ASSERT_EQUALS( p2.getRow(), -1);
......@@ -311,7 +312,6 @@ public:
void test_setQLabFrame2()
{
// Create fictional instrument
const V3D source(0,0,0);
const V3D sample(15, 0, 0);
......@@ -320,24 +320,12 @@ public:
const V3D beam2 = detectorPos - sample;
auto minimalInstrument = ComponentCreationHelper::createMinimalInstrument( source, sample, detectorPos );
// Calculate energy of neutron based on velocity
const double velocity = 1.1 * 10e3; // m/sec
double efixed = 0.5 * Mantid::PhysicalConstants::NeutronMass * velocity * velocity ; // In Joules
efixed = efixed / Mantid::PhysicalConstants::meV;
// Derive distances and angles
const double l1 = beam1.norm();
const double l2 = beam2.norm();
const double scatteringAngle2 = beam2.angle(beam1);
const V3D qLabDir = (beam1/l1) - (beam2/l2);
// Derive the wavelength
std::vector<double> x;
const double microSecsInSec = 1e6;
x.push_back( ( (l1 + l2) / velocity ) * microSecsInSec ); // Make a TOF
std::vector<double> y;
Unit_sptr unitOfLambda = UnitFactory::Instance().create("Wavelength");
unitOfLambda->fromTOF(x, y, l1, l2, scatteringAngle2, 0, efixed, 0);
// Derive QLab for diffraction
const double wavenumber_in_angstrom_times_tof_in_microsec =
......@@ -389,7 +377,7 @@ public:
V3D detPos1 = p1.getDetPos();
// Construct using just Q
Peak p2(inst, Qlab1, detPos1.norm());
Peak p2(inst, Qlab1, boost::optional<double>(detPos1.norm()));
TS_ASSERT( p2.findDetector() );
comparePeaks(p1, p2);
TS_ASSERT_EQUALS( p2.getBankName(), "bank1");
......
......@@ -14,6 +14,8 @@
#include "MantidKernel/V3D.h"
#include "MantidKernel/Strings.h"
#include "MantidKernel/PhysicalConstants.h"
#include "MantidGeometry/Crystal/OrientedLattice.h"
#include "MantidGeometry/Instrument/Goniometer.h"
#include "MantidTestHelpers/NexusTestHelper.h"
#include "MantidTestHelpers/ComponentCreationHelper.h"
#include "MantidAPI/AlgorithmManager.h"
......@@ -360,6 +362,58 @@ public:
TS_ASSERT_EQUALS(coordSystem, pw->getSpecialCoordinateSystem());
}
void test_createPeakHKL()
{
// Create simple fictional instrument
const V3D source(0,0,0);
const V3D sample(15, 0, 0);
const V3D detectorPos(20, 5, 0);
const V3D beam1 = sample - source;
const V3D beam2 = detectorPos - sample;
auto minimalInstrument = ComponentCreationHelper::createMinimalInstrument( source, sample, detectorPos );
// Derive distances and angles
const double l1 = beam1.norm();
const double l2 = beam2.norm();
const V3D qLabDir = (beam1/l1) - (beam2/l2);
const double microSecsInSec = 1e6;
// Derive QLab for diffraction
const double wavenumber_in_angstrom_times_tof_in_microsec =
(Mantid::PhysicalConstants::NeutronMass * (l1 + l2) * 1e-10 * microSecsInSec) /
Mantid::PhysicalConstants::h_bar;
V3D qLab = qLabDir * wavenumber_in_angstrom_times_tof_in_microsec;
Mantid::Geometry::OrientedLattice orientedLattice(1, 1, 1, 90, 90, 90); // U is identity, real and reciprocal lattice vectors are identical.
Mantid::Geometry::Goniometer goniometer; // identity
V3D hkl = qLab / (2 * M_PI); // Given our settings above, this is the simplified relationship between qLab and hkl.
// Now create a peaks workspace around the simple fictional instrument
PeaksWorkspace ws;
ws.setInstrument(minimalInstrument);
ws.mutableSample().setOrientedLattice(&orientedLattice);
ws.mutableRun().setGoniometer(goniometer, false);
// Create the peak
Peak* peak = ws.createPeakHKL(hkl);
/*
Now we check we have made a self - consistent peak
*/
TSM_ASSERT_EQUALS("New peak should have HKL we demanded.", hkl, peak->getHKL());
TSM_ASSERT_EQUALS("New peak should have QLab we expected.", qLab, peak->getQLabFrame());
TSM_ASSERT_EQUALS("QSample and QLab should be identical given the identity goniometer settings.", peak->getQLabFrame(), peak->getQSampleFrame());
auto detector = peak->getDetector();
TSM_ASSERT("No detector", detector);
TSM_ASSERT_EQUALS("This detector id does not match what we expect from the instrument definition", 1, detector->getID());
TSM_ASSERT_EQUALS("Thie detector position is wrong", detectorPos, detector->getPos());
TSM_ASSERT_EQUALS("Goniometer has not been set properly", goniometer.getR(), peak->getGoniometerMatrix());
// Clean up.
delete peak;
}
private:
PeaksWorkspace_sptr createSaveTestPeaksWorkspace()
......
......@@ -9,9 +9,9 @@ namespace {
// logger for the algorithm workspaces
Kernel::Logger g_Log("MDWSTransform");
}
using namespace CnvrtToMD;
/** method to build the Q-coordinates transformation.
*
* @param TargWSDescription -- the class which describes target MD workspace. In
......
......@@ -14,6 +14,11 @@ void setQLabFrame(IPeak &peak, Mantid::Kernel::V3D qLabFrame, double distance) {
// Set the qlab frame
return peak.setQLabFrame(qLabFrame, distance);
}
void setQSampleFrame(IPeak &peak, Mantid::Kernel::V3D qSampleFrame, double distance) {
// Set the qsample frame
return peak.setQSampleFrame(qSampleFrame, distance);
}
}
void export_IPeak()
......@@ -42,8 +47,9 @@ void export_IPeak()
"Using the instrument set in the peak, perform ray tracing to find the exact detector.")
.def("getQSampleFrame", &IPeak::getQSampleFrame, "Return the Q change (of the lattice, k_i - k_f) for this peak."
"The Q is in the Sample frame: the goniometer rotation WAS taken out. ")
.def("setQLabFrame", (void (IPeak::*)(Mantid::Kernel::V3D))&IPeak::setQLabFrame)
.def("setQLabFrame", setQLabFrame, "Set the peak using the peak's position in reciprocal space, in the lab frame.")
.def("setQLabFrame", (void (IPeak::*)(Mantid::Kernel::V3D))&IPeak::setQLabFrame, "Set the peak using the peak's position in reciprocal space, in the lab frame.")
.def("setQLabFrame", setQLabFrame, "Set the peak using the peak's position in reciprocal space, in the lab frame. Detector distance explicitly supplied.") // two argument overload
.def("setQSampleFrame", (void (IPeak::*)(Mantid::Kernel::V3D))&IPeak::setQSampleFrame, "Set the peak using the peak's position in reciprocal space, in the sample frame.")
.def("setQSampleFrame", &IPeak::setQSampleFrame, "Set the peak using the peak's position in reciprocal space, in the sample frame.")
.def("setWavelength", &IPeak::setWavelength, "Set the incident wavelength of the neutron. Calculates the energy from this assuming elastic scattering.")
.def("getWavelength", &IPeak::getWavelength, "Return the incident wavelength")
......
#include "MantidAPI/IPeaksWorkspace.h"
#include "MantidAPI/IPeak.h"
#include "MantidPythonInterface/kernel/Registry/DataItemInterface.h"
#include "MantidPythonInterface/kernel/Converters/PyObjectToV3D.h"
#include <boost/python/class.hpp>
#include <boost/python/return_internal_reference.hpp>
......@@ -8,6 +9,28 @@ using namespace Mantid::API;
using Mantid::PythonInterface::Registry::DataItemInterface;
using namespace boost::python;
namespace {
/// Create a peak via it's HKL value from a list or numpy array
void createPeakHKL(IPeaksWorkspace & self, const object& data)
{
self.createPeakHKL(Mantid::PythonInterface::Converters(data));
}
/// Create a peak via it's QLab value from a list or numpy array
void createPeakQLab(IPeaksWorkspace & self, const object& data)
{
self.createPeak(Mantid::PythonInterface::Converters(data));
}
/// Create a peak via it's QLab value from a list or numpy array
void createPeakQLabWithDistance(IPeaksWorkspace & self, const object& data, double detectorDistance)
{
self.createPeak(Mantid::PythonInterface::Converters(data), distance);
}
}
void export_IPeaksWorkspace()
{
// IPeaksWorkspace class
......@@ -16,7 +39,9 @@ void export_IPeaksWorkspace()
.def("addPeak", &IPeaksWorkspace::addPeak, "Add a peak to the workspace")
.def("removePeak", &IPeaksWorkspace::removePeak, "Remove a peak from the workspace")
.def("getPeak", &IPeaksWorkspace::getPeakPtr, return_internal_reference<>(), "Returns a peak at the given index" )
.def("createPeak", &IPeaksWorkspace::createPeak, return_internal_reference<>(), "Create a Peak and return it")
.def("createPeak", createPeakQLab, return_internal_reference<>(), "Create a Peak and return it from its coordinates in the QLab frame")
.def("createPeak", createPeakQLabWithDistance, return_internal_reference<>(), "Create a Peak and return it from its coordinates in the QLab frame, detector-sample distance explicitly provided")
.def("createPeakHKL", createPeakHKL, return_internal_reference<>(), "Create a Peak and return it from its coordinates in the HKL frame")
.def("hasIntegratedPeaks", &IPeaksWorkspace::hasIntegratedPeaks, "Determine if the peaks have been integrated")
.def("getRun", &IPeaksWorkspace::mutableRun, return_internal_reference<>(),
"Return the Run object for this workspace")
......
......@@ -51,7 +51,7 @@ class IPeaksWorkspaceTest(unittest.TestCase):
# Peaks workspace will not be integrated by default.
self.assertTrue(not pws.hasIntegratedPeaks())
def test_peak_setQLab(self):
def test_peak_setQLabFrame(self):
pws = WorkspaceCreationHelper.createPeaksWorkspace(1)
p = pws.getPeak(0)
try:
......@@ -64,6 +64,19 @@ class IPeaksWorkspaceTest(unittest.TestCase):
except Exception:
self.fail("Tried setQLabFrame with one V3D argument and a double distance")
def test_peak_setQSampleFrame(self):
pws = WorkspaceCreationHelper.createPeaksWorkspace(1)
p = pws.getPeak(0)
try:
p.setQSampleFrame(V3D(1,1,1))
except Exception:
self.fail("Tried setQSampleFrame with one V3D argument")
try:
p.setQSampleFrame(V3D(1,1,1), 1)
except Exception:
self.fail("Tried setQSampleFrame with one V3D argument and a double distance")
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment