Skip to content

Commit

Permalink
Refs #11190. Enabled to set wavelength with HB3A.
Browse files Browse the repository at this point in the history
  • Loading branch information
@wdzhou
wdzhou committed May 19, 2015
1 parent 2dec455 commit 34fe0fd
Show file tree
Hide file tree
Showing 3 changed files with 141 additions and 11 deletions.
5 changes: 5 additions & 0 deletions Code/Mantid/Framework/DataHandling/inc/MantidDataHandling/LoadSpiceXML2DDet.h
View file
Expand Up @@ -102,6 +102,11 @@ class DLLExport LoadSpiceXML2DDet : public API::Algorithm {
/// Load instrument
void loadInstrument(API::MatrixWorkspace_sptr matrixws,
const std::string &idffilename);

/// Get wavelength from workspace
bool getHB3AWavelength(API::MatrixWorkspace_sptr dataws, double &wavelength);

void setXtoLabQ(API::MatrixWorkspace_sptr dataws, const double &wavelength);
};

} // namespace DataHandling
Expand Down
100 changes: 95 additions & 5 deletions Code/Mantid/Framework/DataHandling/src/LoadSpiceXML2DDet.cpp
View file
Expand Up @@ -242,6 +242,15 @@ void LoadSpiceXML2DDet::exec() {
std::vector<SpiceXMLNode> vec_xmlnode;
parseSpiceXML(xmlfilename, vec_xmlnode);

// Spice data in table workspace
std::string spicetablewsname = getPropertyValue("SpiceTableWorkspace");
ITableWorkspace_sptr spicetablews;
int ptnumber(0);
if (spicetablewsname.size() > 0) {
spicetablews = getProperty("SpiceTableWorkspace");
ptnumber = getProperty("PtNumber");
}

// Create output workspace
MatrixWorkspace_sptr outws;
if (loadinstrument) {
Expand All @@ -252,15 +261,21 @@ void LoadSpiceXML2DDet::exec() {
false);
}

std::string spicetablewsname = getPropertyValue("SpiceTableWorkspace");
// Set up log
if (spicetablewsname.size() > 0) {
ITableWorkspace_sptr spicetablews = getProperty("SpiceTableWorkspace");
int ptnumber = getProperty("PtNumber");
setupSampleLogFromSpiceTable(outws, spicetablews, ptnumber);
}

if (loadinstrument)
if (loadinstrument) {
loadInstrument(outws, idffilename);
if (spicetablewsname.size() > 0) {
double wavelength;
bool has_wavelength = getHB3AWavelength(outws, wavelength);
if (has_wavelength) {
setXtoLabQ(outws, wavelength);
}
}
}

setProperty("OutputWorkspace", outws);
}
Expand Down Expand Up @@ -383,7 +398,7 @@ MatrixWorkspace_sptr LoadSpiceXML2DDet::createMatrixWorkspace(
if (loadinstrument) {
size_t numspec = numpixelx * numpixely;
outws = boost::dynamic_pointer_cast<MatrixWorkspace>(
WorkspaceFactory::Instance().create("Workspace2D", numspec, 1, 1));
WorkspaceFactory::Instance().create("Workspace2D", numspec, 2, 1));
} else {
outws = boost::dynamic_pointer_cast<MatrixWorkspace>(
WorkspaceFactory::Instance().create("Workspace2D", numpixely, numpixelx,
Expand Down Expand Up @@ -547,6 +562,81 @@ void LoadSpiceXML2DDet::setupSampleLogFromSpiceTable(
return;
}

bool LoadSpiceXML2DDet::getHB3AWavelength(MatrixWorkspace_sptr dataws,
double &wavelength) {
bool haswavelength(false);
wavelength = -1;

if (dataws->run().hasProperty("_m1")) {
Kernel::TimeSeriesProperty<double> *ts =
dynamic_cast<Kernel::TimeSeriesProperty<double> *>(
dataws->run().getProperty("_m1"));
if (ts && ts->size() > 0) {
double m1pos = ts->valuesAsVector()[0];
if (fabs(m1pos - (-25.870000)) < 0.2) {
wavelength = 1.003;
haswavelength = true;
} else {
g_log.warning() << "m1 position " << m1pos
<< " does not have defined mapping to "
<< "wavelength."
<< "\n";
}
} else if (!ts) {
g_log.warning("Log _m1 is not TimeSeriesProperty. Treat it as a single "
"value property.");
double m1pos = atof(dataws->run().getProperty("_m1")->value().c_str());
if (fabs(m1pos - (-25.870000)) < 0.2) {
wavelength = 1.003;
haswavelength = true;
} else {
g_log.warning() << "m1 position " << m1pos
<< " does not have defined mapping to "
<< "wavelength."
<< "\n";
}
} else {
g_log.error("Log _m1 is empty.");
}
} else {
g_log.warning() << "No _m1 log is found."
<< "\n";
}

if (!haswavelength)
g_log.warning("No wavelength is setup!");
else
g_log.notice() << "[DB] Wavelength = " << wavelength << "\n";

return haswavelength;
}

void LoadSpiceXML2DDet::setXtoLabQ(API::MatrixWorkspace_sptr dataws,
const double &wavelength) {
// Geometry information
// Kernel::V3D sourcepos = dataws->getInstrument()->getSource()->getPos();
// Kernel::V3D samplepos = dataws->getInstrument()->getSample()->getPos();
// Kernel::V3D sample_source_vector = sourcepos - samplepos;

size_t numspec = dataws->getNumberHistograms();
for (size_t iws = 0; iws < numspec; ++iws) {
/*
Kernel::V3D detpos = dataws->getDetector(iws)->getPos();
Kernel::V3D detpos_sample_vector = detpos - samplepos;
double scattering_rad = detpos_sample_vector.angle(sample_source_vector);
*/

// double ki = 4*sin(scattering_rad*0.5)*M_PI/wavelength;
double ki = 2. * M_PI / wavelength;
dataws->dataX(iws)[0] = ki;
dataws->dataX(iws)[1] = ki + 0.00001;
}

dataws->getAxis(0)->setUnit("Momentum");

return;
}

//----------------------------------------------------------------------------------------------
/** Load instrument
* @brief LoadSpiceXML2DDet::loadInstrument
Expand Down
47 changes: 41 additions & 6 deletions Code/Mantid/Framework/DataHandling/test/LoadSpiceXML2DDetTest.h
View file
Expand Up @@ -120,14 +120,15 @@ class LoadSpiceXML2DDetTest : public CxxTest::TestSuite {
boost::make_shared<DataObjects::TableWorkspace>());
datatablews->addColumn("int", "Pt.");
datatablews->addColumn("double", "2theta");
datatablews->addColumn("double", "m1");
AnalysisDataService::Instance().addOrReplace("SpiceDataTable", datatablews);

TableRow row0 = datatablews->appendRow();
row0 << 1 << 15.0;
row0 << 1 << 15.0 << -25.8;
TableRow row1 = datatablews->appendRow();
row1 << 2 << -15.0;
row1 << 2 << -15.0 << -25.8;
TableRow row2 = datatablews->appendRow();
row2 << 3 << 0.0;
row2 << 3 << 0.0 << -25.8;

// Test 2theta = 15 degree
LoadSpiceXML2DDet loader;
Expand Down Expand Up @@ -169,6 +170,16 @@ class LoadSpiceXML2DDetTest : public CxxTest::TestSuite {
Kernel::V3D detlast = outws->getDetector(256 * 256 - 1)->getPos();
Kernel::V3D detmiddle =
outws->getDetector(256 / 2 * 256 + 256 / 2)->getPos();
Kernel::V3D detedgemiddle = outws->getDetector(256 * 256 / 2)->getPos();
std::cout << "\n";
std::cout << "Sample position: " << sample.X() << ", " << sample.Y() << ", "
<< sample.Z() << "\n";
std::cout << "Source position: " << source.X() << ", " << source.Y() << ", "
<< source.Z() << "\n";
std::cout << "Det Edge Middle: " << detedgemiddle.X() << ", "
<< detedgemiddle.Y() << ", " << detedgemiddle.Z() << "\n";
std::cout << "Det (1, 1) : " << det0pos.X() << ", " << det0pos.Y()
<< ", " << det0pos.Z() << "\n";

// center of the detector must be negative
TS_ASSERT_LESS_THAN(detmiddle.X(), 0.0);
Expand All @@ -187,10 +198,11 @@ class LoadSpiceXML2DDetTest : public CxxTest::TestSuite {

// 2theta value
Kernel::V3D sample_source = sample - source;
std::cout << "Sample - Source = " << sample_source.X() << ", "
<< sample_source.Y() << ", " << sample_source.Z() << "\n";

Kernel::V3D det0_sample = det0pos - sample;
double twotheta0 = det0_sample.angle(sample_source) * 180. / 3.14159265;
// FIXME - Verify this with Huibo!
TS_ASSERT_DELTA(twotheta0, 11.6252, 0.0001);

Kernel::V3D detTL_sample = detlast0 - sample;
Expand All @@ -199,6 +211,7 @@ class LoadSpiceXML2DDetTest : public CxxTest::TestSuite {

Kernel::V3D det255_sample = det255pos - sample;
double twotheta255 = det255_sample.angle(sample_source) * 180. / 3.14159265;
TS_ASSERT_DELTA(twotheta255, 19.5328, 0.0001);

Kernel::V3D detlast_sample = detlast - sample;
double twothetalast =
Expand All @@ -211,10 +224,17 @@ class LoadSpiceXML2DDetTest : public CxxTest::TestSuite {
detmid_sample.angle(sample_source) * 180. / 3.1415926;
TS_ASSERT_DELTA(twotheta_middle, 15.0, 0.02);

Kernel::V3D detedgemid_sample = detedgemiddle - sample;
double twotheta_edgemiddle =
detedgemid_sample.angle(sample_source) * 180. / 3.14159265;
TS_ASSERT_DELTA(twotheta_edgemiddle, 10.8864, 0.0001);

TS_ASSERT_LESS_THAN(twotheta0, twotheta_middle);
TS_ASSERT_LESS_THAN(twotheta_middle, twothetalast);

//-------------------------------------------------------------------------------
// Case: 2theta = 0
//-------------------------------------------------------------------------------
LoadSpiceXML2DDet loader2;
loader2.initialize();

Expand Down Expand Up @@ -243,6 +263,8 @@ class LoadSpiceXML2DDetTest : public CxxTest::TestSuite {
Kernel::V3D detlast2 = outws2->getDetector(256 * 256 - 1)->getPos();
Kernel::V3D detmiddle2 =
outws2->getDetector(256 / 2 * 256 + 256 / 2)->getPos();
std::cout << "Case 2: Middle Det " << detmiddle2.X() << ", "
<< detmiddle2.Y() << ", " << detmiddle2.Z() << "\n";

// detector distance
double dist0b = sample2.distance(det0pos2);
Expand All @@ -261,21 +283,34 @@ class LoadSpiceXML2DDetTest : public CxxTest::TestSuite {

detmid_sample = detmiddle2 - sample2;
twotheta_middle = detmid_sample.angle(sample_source) * 180. / 3.1415926;
// FIXME - Justify with Huibo
TS_ASSERT_DELTA(twotheta_middle, 0.00, 0.03);
TS_ASSERT_DELTA(twotheta_middle, 0.0228, 0.0001);

det0_sample = det0pos2 - sample2;
double twotheta0_2 = det0_sample.angle(sample_source) * 180. / 3.14159265;
std::cout << "Case 2: 2theta at (1, 1) " << twotheta0_2 << "\n";

detTL_sample = detlast0_2 - sample2;
twotheta_tl = detTL_sample.angle(sample_source) * 180. / 3.14159265;

det255_sample = det255pos2 - sample2;
double twotheta255_2 =
det255_sample.angle(sample_source) * 180. / 3.14159265;
std::cout << "Case 2: 2theta at (1, 256) " << twotheta255_2 << "\n";

detlast_sample = detlast2 - sample2;
twothetalast = detlast_sample.angle(sample_source) * 180. / 3.14159265;
std::cout << "Case 2: 2theta at (256, 256) " << twothetalast << "\n";

Kernel::V3D det_edgemiddel_left =
outws2->getDetector(256 * 256 / 2)->getPos();
Kernel::V3D det_edgemiddel_right =
outws2->getDetector(256 * 256 / 2 + 255)->getPos();
Kernel::V3D det_edgemiddle_left_sample = det_edgemiddel_left - sample;
Kernel::V3D det_edgemiddle_right_sample = det_edgemiddel_right - sample;
double inplane1 = det_edgemiddle_left_sample.angle(sample_source);
double inplane2 = det_edgemiddle_right_sample.angle(sample_source);
std::cout << "In plain left side = " << inplane1 * 180 / 3.1415926 << "\n";
std::cout << "In plain right side = " << inplane2 * 180 / 3.1415926 << "\n";

TS_ASSERT_DELTA(twotheta0_2, twothetalast, 0.00001);
TS_ASSERT_DELTA(twotheta0_2, twotheta_tl, 0.00001);
Expand Down

0 comments on commit 34fe0fd

Please sign in to comment.