diff --git a/Code/Mantid/Framework/API/src/ExperimentInfo.cpp b/Code/Mantid/Framework/API/src/ExperimentInfo.cpp index d7a43ffc2ec1..5f7b9d32f4c4 100644 --- a/Code/Mantid/Framework/API/src/ExperimentInfo.cpp +++ b/Code/Mantid/Framework/API/src/ExperimentInfo.cpp @@ -141,6 +141,7 @@ const std::string ExperimentInfo::toString() const { */ void ExperimentInfo::setInstrument(const Instrument_const_sptr &instr) { if (instr->isParametrized()) { + g_log.notice("[DB] ExperimentInfo.setInstrument as being isParametrized."); sptr_instrument = instr->baseInstrument(); m_parmap = instr->getParameterMap(); } else { diff --git a/Code/Mantid/Framework/DataHandling/src/LoadSpiceXML2DDet.cpp b/Code/Mantid/Framework/DataHandling/src/LoadSpiceXML2DDet.cpp index 40c1e500b241..035a85e1e63b 100644 --- a/Code/Mantid/Framework/DataHandling/src/LoadSpiceXML2DDet.cpp +++ b/Code/Mantid/Framework/DataHandling/src/LoadSpiceXML2DDet.cpp @@ -265,6 +265,24 @@ void LoadSpiceXML2DDet::exec() { if (spicetablewsname.size() > 0) { setupSampleLogFromSpiceTable(outws, spicetablews, ptnumber); } + else + { + // Set up 2theta from _2theta + if (outws->run().hasProperty("_2theta")) + { + Kernel::DateAndTime anytime(1000); + double logvalue = atof(outws->run().getProperty("_2theta")->value().c_str()); + g_log.information() << "Set 2theta from _2thea with value " << logvalue << "\n"; + TimeSeriesProperty *newlogproperty = + new TimeSeriesProperty("2theta"); + newlogproperty->addValue(anytime, logvalue); + outws->mutableRun().addProperty(newlogproperty); + } + else + { + g_log.warning("No 2theta is set up for loading instrument."); + } + } if (loadinstrument) { loadInstrument(outws, idffilename); @@ -613,20 +631,9 @@ bool LoadSpiceXML2DDet::getHB3AWavelength(MatrixWorkspace_sptr dataws, 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; diff --git a/Code/Mantid/Framework/DataHandling/test/LoadSpiceXML2DDetTest.h b/Code/Mantid/Framework/DataHandling/test/LoadSpiceXML2DDetTest.h index 9f6a8781516d..055090fee295 100644 --- a/Code/Mantid/Framework/DataHandling/test/LoadSpiceXML2DDetTest.h +++ b/Code/Mantid/Framework/DataHandling/test/LoadSpiceXML2DDetTest.h @@ -359,6 +359,86 @@ class LoadSpiceXML2DDetTest : public CxxTest::TestSuite { AnalysisDataService::Instance().remove("Exp0335_S0038B"); AnalysisDataService::Instance().remove("Exp0335_S0038"); } + + + void test_LoadHB3AXMLInstrumentNoTable() + { + // Test 2theta = 15 degree + LoadSpiceXML2DDet loader; + loader.initialize(); + + const std::string filename("HB3A_exp355_scan0001_0522.xml"); + TS_ASSERT_THROWS_NOTHING(loader.setProperty("Filename", filename)); + TS_ASSERT_THROWS_NOTHING( + loader.setProperty("OutputWorkspace", "Exp0335_S0038")); + std::vector sizelist(2); + sizelist[0] = 256; + sizelist[1] = 256; + loader.setProperty("DetectorGeometry", sizelist); + loader.setProperty("LoadInstrument", true); + + loader.execute(); + TS_ASSERT(loader.isExecuted()); + + // Get data + MatrixWorkspace_sptr outws = boost::dynamic_pointer_cast( + AnalysisDataService::Instance().retrieve("Exp0335_S0038")); + TS_ASSERT(outws); + TS_ASSERT_EQUALS(outws->getNumberHistograms(), 256 * 256); + + // Value + TS_ASSERT_DELTA(outws->readY(255)[0], 1.0, 0.0001); + TS_ASSERT_DELTA(outws->readY(253 * 256 + 9)[0], 1.0, 0.00001); + + // Instrument + TS_ASSERT(outws->getInstrument()); + + Kernel::V3D source = outws->getInstrument()->getSource()->getPos(); + Kernel::V3D sample = outws->getInstrument()->getSample()->getPos(); + Kernel::V3D det0pos = outws->getDetector(0)->getPos(); + Kernel::V3D det255pos = outws->getDetector(255)->getPos(); + Kernel::V3D detlast0 = outws->getDetector(256 * 255)->getPos(); + 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); + + // detector distance + double dist0 = sample.distance(det0pos); + double dist255 = sample.distance(det255pos); + double distlast = sample.distance(detlast); + double distlast0 = sample.distance(detlast0); + double distmiddle = sample.distance(detmiddle); + + TS_ASSERT_DELTA(dist0, dist255, 0.0001); + TS_ASSERT_DELTA(dist0, distlast, 0.0001); + TS_ASSERT_DELTA(dist0, distlast0, 0.0001); + TS_ASSERT_DELTA(distmiddle, 0.3518, 0.000001); + + // 2theta value + Kernel::V3D sample_source = sample - source; + std::cout << "Sample - Source = " << sample_source.X() << ", " + << sample_source.Y() << ", " << sample_source.Z() << "\n"; + + + Kernel::V3D detmid_sample = detmiddle - sample; + double twotheta_middle = + detmid_sample.angle(sample_source) * 180. / 3.1415926; + TS_ASSERT_DELTA(twotheta_middle, 42.70975, 0.02); + } + }; #endif /* MANTID_DATAHANDLING_LOADSPICEXML2DDETTEST_H_ */ diff --git a/Code/Mantid/Framework/MDAlgorithms/inc/MantidMDAlgorithms/ConvertCWSDExpToMomentum.h b/Code/Mantid/Framework/MDAlgorithms/inc/MantidMDAlgorithms/ConvertCWSDExpToMomentum.h index 560f2ea737e4..914cdbd278c1 100644 --- a/Code/Mantid/Framework/MDAlgorithms/inc/MantidMDAlgorithms/ConvertCWSDExpToMomentum.h +++ b/Code/Mantid/Framework/MDAlgorithms/inc/MantidMDAlgorithms/ConvertCWSDExpToMomentum.h @@ -33,11 +33,29 @@ namespace MDAlgorithms { File change history is stored at: Code Documentation is available at: */ -class DLLExport ConvertCWSDExpToMomentum : API::Algorithm { +class DLLExport ConvertCWSDExpToMomentum : public API::Algorithm { public: ConvertCWSDExpToMomentum(); virtual ~ConvertCWSDExpToMomentum(); + /// Algorithm's name + virtual const std::string name() const { + return "ConvertCWSDExpToMomentum"; + } + + /// Summary of algorithms purpose + virtual const std::string summary() const { + return "Load and convert a set of files in an HB3A experiment."; + } + + /// Algorithm's version + virtual int version() const { return (1); } + + /// Algorithm's category for identification + virtual const std::string category() const { + return "Diffraction;DataHandling\\Text"; + } + private: void init(); void exec(); @@ -48,6 +66,9 @@ class DLLExport ConvertCWSDExpToMomentum : API::Algorithm { const detid_t &startdetid, const int runnumber); + void convertToMomentum(const std::vector &detPos, const double &wavelength, + std::vector &qSample); + API::IMDEventWorkspace_sptr createExperimentMDWorkspace(); bool getInputs(std::string &errmsg); @@ -61,6 +82,7 @@ class DLLExport ConvertCWSDExpToMomentum : API::Algorithm { API::ITableWorkspace_sptr m_expDataTableWS; API::ITableWorkspace_sptr m_detectorListTableWS; API::IMDEventWorkspace_sptr m_outputWS; + Geometry::Instrument_sptr m_virtualInstrument; Kernel::V3D m_samplePos; Kernel::V3D m_sourcePos; @@ -68,10 +90,9 @@ class DLLExport ConvertCWSDExpToMomentum : API::Algorithm { size_t m_iColFilename; size_t m_iColStartDetID; - size_t m_numBins; - std::vector m_extentMins; std::vector m_extentMaxs; + std::vector m_numBins; }; } // namespace MDAlgorithms diff --git a/Code/Mantid/Framework/MDAlgorithms/src/ConvertCWSDExpToMomentum.cpp b/Code/Mantid/Framework/MDAlgorithms/src/ConvertCWSDExpToMomentum.cpp index d7e1dd42be08..80babb5514a8 100644 --- a/Code/Mantid/Framework/MDAlgorithms/src/ConvertCWSDExpToMomentum.cpp +++ b/Code/Mantid/Framework/MDAlgorithms/src/ConvertCWSDExpToMomentum.cpp @@ -3,6 +3,7 @@ #include "MantidKernel/ArrayProperty.h" #include "MantidGeometry/Instrument.h" #include "MantidDataObjects/MDEventFactory.h" +#include "MantidAPI/ExperimentInfo.h" #include "MantidGeometry/Instrument/ComponentHelper.h" using namespace Mantid::API; @@ -13,6 +14,8 @@ using namespace Mantid::Geometry; namespace Mantid { namespace MDAlgorithms { + DECLARE_ALGORITHM(ConvertCWSDExpToMomentum) + //---------------------------------------------------------------------------------------------- /** Constructor */ @@ -36,6 +39,10 @@ void ConvertCWSDExpToMomentum::init() { "DetectorTableWorkspace", "", Direction::Input), "Name of table workspace containing all the detectors."); + declareProperty(new WorkspaceProperty( + "OutputWorkspace", "", Direction::Output), + "Name of MDEventWorkspace containing all experimental data."); + declareProperty(new ArrayProperty("SourcePosition"), "A vector of 3 doubles for position of source."); @@ -44,6 +51,7 @@ void ConvertCWSDExpToMomentum::init() { declareProperty(new ArrayProperty("PixelDimension"), "A vector of 8 doubles to determine a cubic pixel's size."); + } /** @@ -56,13 +64,17 @@ void ConvertCWSDExpToMomentum::exec() { if (!inputvalid) throw std::runtime_error(errmsg); + g_log.notice("[DB] Inputs are examined."); + // Create output MDEventWorkspace - std::vector vecDetPos; - std::vector vecDetID; m_outputWS = createExperimentMDWorkspace(); // Add MDEventWorkspace addMDEvents(); + + setProperty("OutputWorkspace", m_outputWS); + + return; } API::IMDEventWorkspace_sptr @@ -73,10 +85,20 @@ ConvertCWSDExpToMomentum::createExperimentMDWorkspace() { parseDetectorTable(vec_detpos, vec_detid); // FIXME - Should set create a solid instrument - Geometry::Instrument_sptr virtualinstrument = Geometry::ComponentHelper::createVirtualInstrument(m_sourcePos, - m_samplePos, - vec_detpos, - vec_detid); + m_virtualInstrument = Geometry::ComponentHelper::createVirtualInstrument(m_sourcePos, m_samplePos, + vec_detpos, + vec_detid); + if (!m_virtualInstrument) + throw std::runtime_error("Failed to create virtual instrument."); + g_log.notice() << "[DB] Virtual Instrument has " << m_virtualInstrument->getDetectorIDs().size() + << "Detectors\n"; + + + // Create workspace + // FIXME - Should I give out a better value??? + m_extentMins.resize(3, -10.0); + m_extentMaxs.resize(3, 10.0); + m_numBins.resize(3, 100); size_t m_nDimensions = 3; IMDEventWorkspace_sptr mdws = @@ -95,17 +117,22 @@ ConvertCWSDExpToMomentum::createExperimentMDWorkspace() { dimensionNames[1] = "Q_sample_y"; dimensionNames[2] = "Q_sample_z"; Mantid::Kernel::SpecialCoordinateSystem coordinateSystem = Mantid::Kernel::QSample; + // Add dimensions + if (m_extentMins.size() != 3 || m_extentMaxs.size() != 3 || m_numBins.size() != 3) + throw std::runtime_error("The range and number of bins are not set up correctly to create MDEventWorkspace."); + else + for (size_t d = 0; d < 3; ++d) + g_log.debug() << "Direction " << d << ", Range = " << m_extentMins[d] + << ", " << m_extentMaxs[d] << "\n"; + for (size_t i = 0; i < m_nDimensions; ++i) { std::string id = vec_ID[i]; - std::string name = vec_name[i]; + std::string name = dimensionNames[i]; // FIXME - Unit of momentum? std::string units = "???"; // int nbins = 100; - for (size_t d = 0; d < 3; ++d) - g_log.debug() << "Direction " << d << ", Range = " << m_extentMins[d] - << ", " << m_extentMaxs[d] << "\n"; mdws->addDimension( Geometry::MDHistoDimension_sptr(new Geometry::MDHistoDimension( id, name, units, static_cast(m_extentMins[i]), @@ -158,20 +185,35 @@ void ConvertCWSDExpToMomentum::addMDEvents() { void ConvertCWSDExpToMomentum::convertSpiceMatrixToMomentumMDEvents( MatrixWorkspace_sptr dataws, const detid_t &startdetid, const int runnumber) { + + // Creates a new instance of the MDEventInserter. + MDEventWorkspace, 3>::sptr mdws_mdevt_3 = + boost::dynamic_pointer_cast, 3> >(m_outputWS); + MDEventInserter, 3>::sptr> inserter( + mdws_mdevt_3); + size_t numspec = dataws->getNumberHistograms(); for (size_t iws = 0; iws < numspec; ++iws) { + // Get detector positions and signal double signal = dataws->readY(iws)[0]; double error = dataws->readE(iws)[0]; double wavelength = dataws->readX(iws)[0]; - + // Create the MDEvent Kernel::V3D detpos = dataws->getDetector(iws)->getPos(); - std::vector momentum = - convertToMomentum(detpos, wavelength, goniometersetup); + std::vector momentum(3); + convertToMomentum(detpos, wavelength, momentum); detid_t detid = dataws->getDetector(iws)->getID() + startdetid; - - // Create the MDEvent - m_outputWS->addExperimentInfo(); + // Insert + inserter.insertMDEvent(static_cast(signal), static_cast(error*error), + static_cast(runnumber), detid, + momentum.data()); + // m_outputWS->addExperimentInfo(); } + + ExperimentInfo_sptr expinfo = boost::make_shared(); + expinfo->setInstrument(m_virtualInstrument); + m_outputWS->addExperimentInfo(expinfo); + } /// Examine input @@ -231,6 +273,36 @@ bool ConvertCWSDExpToMomentum::getInputs(std::string &errmsg) { return (errmsg.size() > 0); } +void ConvertCWSDExpToMomentum::convertToMomentum(const std::vector &detPos, const double &wavelength, std::vector &qSample) +{ + // TODO/FIXME : Implement ASAP! + + return; +} + +API::MatrixWorkspace_sptr ConvertCWSDExpToMomentum::loadSpiceData(const std::string &filename, bool &loaded, std::string &errmsg) +{ + IAlgorithm_sptr loader = createChildAlgorithm("LoadSpiceXML2DDet"); + loader->initialize(); + loader->setProperty("Filename", filename); + std::vector sizelist(2); + sizelist[0] = 256; + sizelist[1] = 256; + loader->setProperty("DetectorGeometry", sizelist); + loader->setProperty("LoadInstrument", true); + + loader->execute(); + + API::MatrixWorkspace_sptr dataws = loader->getProperty("OutputWorkspace"); + if (dataws) + loaded = true; + else + loaded = false; + + return dataws; +} + + void ConvertCWSDExpToMomentum::parseDetectorTable( std::vector &vec_detpos, std::vector &vec_detid) { // Set vectors' sizes diff --git a/Code/Mantid/Framework/MDAlgorithms/test/ConvertCWSDExpToMomentumTest.h b/Code/Mantid/Framework/MDAlgorithms/test/ConvertCWSDExpToMomentumTest.h index 0162d19e43f1..b532177114a5 100644 --- a/Code/Mantid/Framework/MDAlgorithms/test/ConvertCWSDExpToMomentumTest.h +++ b/Code/Mantid/Framework/MDAlgorithms/test/ConvertCWSDExpToMomentumTest.h @@ -4,9 +4,16 @@ #include #include "MantidMDAlgorithms/ConvertCWSDExpToMomentum.h" +#include "MantidAPI/ITableWorkspace.h" +#include "MantidAPI/TableRow.h" +#include "MantidDataObjects/TableWorkspace.h" +#include "MantidAPI/IMDIterator.h" using Mantid::MDAlgorithms::ConvertCWSDExpToMomentum; +using namespace Mantid; using namespace Mantid::API; +using namespace Mantid::DataObjects; +using namespace Mantid::Kernel; class ConvertCWSDExpToMomentumTest : public CxxTest::TestSuite { public: @@ -16,13 +23,96 @@ class ConvertCWSDExpToMomentumTest : public CxxTest::TestSuite { static void destroySuite( ConvertCWSDExpToMomentumTest *suite ) { delete suite; } - void test_Something() + void test_Init() { - TSM_ASSERT( "You forgot to write a test!", 0); + // Create tableworkspaces as inputs + generateTestInputs(); + + ConvertCWSDExpToMomentum testalg; + testalg.initialize(); + TS_ASSERT(testalg.isInitialized()); + } + + void test_LoadConvert1File() + { + // Init and set up + ConvertCWSDExpToMomentum testalg; + testalg.initialize(); + + testalg.setProperty("InputWorkspace", "DataFileTable"); + testalg.setProperty("DetectorTableWorkspace", "DetectorTable"); + testalg.setProperty("SourcePosition", m_sourcePos); + testalg.setProperty("SamplePosition", m_samplePos); + testalg.setProperty("PixelDimension", m_pixelDimension); + testalg.setProperty("OutputWorkspace", "QSampleMDEvents"); + + testalg.execute(); + TS_ASSERT(testalg.isExecuted()); + + API::IMDEventWorkspace_sptr outws = boost::dynamic_pointer_cast( + AnalysisDataService::Instance().retrieve("QSampleMDEvents")); + TS_ASSERT(outws); + + IMDIterator *mditer = outws->createIterator(); + TS_ASSERT_EQUALS(mditer->getNumEvents(), 256*256); + + } + + + +private: + API::ITableWorkspace_sptr m_dataTableWS; + API::ITableWorkspace_sptr m_detectorTableWS; + std::vector m_sourcePos; + std::vector m_samplePos; + std::vector m_pixelDimension; + + void generateTestInputs() + { + // Create data table + DataObjects::TableWorkspace_sptr datatable = boost::make_shared(); + datatable->addColumn("int", "Scan No"); + datatable->addColumn("int", "Pt. No"); + datatable->addColumn("str", "File Name"); + datatable->addColumn("int", "Starting DetID"); + TableRow row0 = datatable->appendRow(); + row0 << 1 << 522 << "HB3A_exp355_scan0001_0522.xml" << 256*256; + m_dataTableWS = boost::dynamic_pointer_cast(datatable); + TS_ASSERT(m_dataTableWS); + + + // Create detector table + DataObjects::TableWorkspace_sptr dettable = boost::make_shared(); + dettable->addColumn("int", "DetID"); + dettable->addColumn("double", "X"); + dettable->addColumn("double", "Y"); + dettable->addColumn("double", "Z"); + for (size_t i = 0; i < 256; ++i) + { + TableRow detrow = dettable->appendRow(); + double x = 0.38+static_cast(i-128)*0.001; + double y = 0; + double z = 0.38+static_cast(i-128)*0.001; + detrow << static_cast(i) + 1 << x << y << z; + } + + m_detectorTableWS = boost::dynamic_pointer_cast(dettable); + TS_ASSERT(m_detectorTableWS); + + AnalysisDataService::Instance().addOrReplace("DataFileTable", m_dataTableWS); + AnalysisDataService::Instance().addOrReplace("DetectorTable", m_detectorTableWS); + + // Source and sample position + m_sourcePos.resize(3, 0.0); + m_sourcePos[2] = -2; + + m_samplePos.resize(3, 0.0); + + m_pixelDimension.resize(8, 0.0); } }; -#endif /* MANTID_MDALGORITHMS_CONVERTCWSDEXPTOMOMENTUMTEST_H_ */ \ No newline at end of file +#endif /* MANTID_MDALGORITHMS_CONVERTCWSDEXPTOMOMENTUMTEST_H_ */