diff --git a/Code/Mantid/Framework/Algorithms/src/EditInstrumentGeometry.cpp b/Code/Mantid/Framework/Algorithms/src/EditInstrumentGeometry.cpp index 7a9310403471a08ea28ade16d591ed4f14d85d1f..0e5aa604216768e39fc994b89d0157f40be45716 100644 --- a/Code/Mantid/Framework/Algorithms/src/EditInstrumentGeometry.cpp +++ b/Code/Mantid/Framework/Algorithms/src/EditInstrumentGeometry.cpp @@ -26,6 +26,7 @@ There are some limitations of this algorithm. #include "MantidAlgorithms/EditInstrumentGeometry.h" #include "MantidAPI/IAlgorithm.h" #include "MantidGeometry/IDetector.h" +#include "MantidGeometry/IObjComponent.h" #include "MantidAPI/ISpectrum.h" #include "MantidKernel/System.h" #include "MantidKernel/ArrayProperty.h" @@ -89,9 +90,6 @@ namespace Algorithms declareProperty(new WorkspaceProperty<>("Workspace", "", Direction::InOut), "Workspace to edit the detector information"); - // Instrument Name - declareProperty("InstrumentName", "GenericPowder", "Name of the instrument. "); - // L1 declareProperty("PrimaryFlightPath", EMPTY_DBL(), "Primary flight path L1 of the powder diffractomer. "); @@ -120,8 +118,11 @@ namespace Algorithms declareProperty(new ArrayProperty<int>("DetectorIDs"), "User specified detector IDs of the spectra. " "Number of specified detector IDs must be either zero or number of histogram"); - // Indicator if it is a new instrument - declareProperty("NewInstrument", false, "Add detectors information from scratch"); + // Instrument Name + declareProperty("InstrumentName", "", "Name of the newly built instrument. If left empty, " + "the original instrument will be used. "); + + return; } template <typename NumT> @@ -188,24 +189,30 @@ namespace Algorithms */ void EditInstrumentGeometry::exec() { - // lots of things have to do with the input workspace + // Lots of things have to do with the input workspace MatrixWorkspace_sptr workspace = getProperty("Workspace"); Geometry::Instrument_const_sptr originstrument = workspace->getInstrument(); - const bool newinstrument = getProperty("NewInstrument"); - // get the primary flight path + // Get and check the primary flight path double l1 = this->getProperty("PrimaryFlightPath"); - // Check validity of the input properties if (isEmpty(l1)) { - Kernel::V3D sourcepos = originstrument->getSource()->getPos(); - l1 = fabs(sourcepos.Z()); - g_log.information() << "Retrieve L1 from input data workspace. Source position = " << sourcepos.X() - << ", " << sourcepos.Y() << ", " << sourcepos.Z() << std::endl; + // Use the original L1 + if (!originstrument) + { + std::string errmsg("It is not supported that L1 is not given, ", + "while there is no instrument associated to input workspace."); + g_log.error(errmsg); + throw std::runtime_error(errmsg); + } + Geometry::IObjComponent_const_sptr source = originstrument->getSource(); + Geometry::IObjComponent_const_sptr sample = originstrument->getSample(); + l1 = source->getDistance(*sample); + g_log.information() << "Retrieve L1 from input data workspace. \n"; } g_log.information() << "Using L1 = " << l1 << "\n"; - // get spectra number in case they are in a funny order + // Get spectra number in case they are in a funny order std::vector<int32_t> specids = this->getProperty("SpectrumIDs"); if (specids.empty()) // they are using the order of the input workspace { @@ -216,23 +223,32 @@ namespace Algorithms } } - // get the detector ids - empsy means ignore it + // Get the detector ids - empsy means ignore it const vector<int> vec_detids = getProperty("DetectorIDs"); const bool renameDetID(!vec_detids.empty()); - // individual detector geometries + // Get individual detector geometries ordered by input spectrum IDs const std::vector<double> l2s = this->getProperty("L2"); const std::vector<double> tths = this->getProperty("Polar"); std::vector<double> phis = this->getProperty("Azimuthal"); + // empty list of L2 and 2-theta is not allowed + if (l2s.empty()) + { + throw std::runtime_error("User must specify L2 for all spectra. "); + } + if (tths.empty()) + { + throw std::runtime_error("User must specify 2theta for all spectra."); + } + // empty list of phi means that they are all zero if (phis.empty()) { phis.assign(l2s.size(), 0.); } - // only did work down to here - + // Validate for (size_t ib = 0; ib < l2s.size(); ib ++) { g_log.information() << "Detector " << specids[ib] << " L2 = " << l2s[ib] << " 2Theta = " << tths[ib] << std::endl; @@ -252,28 +268,18 @@ namespace Algorithms // Keep original instrument and set the new instrument, if necessary spec2index_map *spec2indexmap = workspace->getSpectrumToWorkspaceIndexMap(); - // ??? Condition: spectrum has 1 and only 1 detector + // ??? Condition: spectrum has 1 and only 1 detector size_t nspec = workspace->getNumberHistograms(); std::vector<bool> wsindexsetflag(nspec, false); + + // Initialize another set of L2/2-theta/Phi/DetectorIDs vector ordered by workspace index std::vector<bool> storable(nspec, false); std::vector<double> storL2s(nspec, 0.); std::vector<double> stor2Thetas(nspec, 0.); std::vector<double> storPhis(nspec, 0.); vector<int> storDetIDs(nspec, 0); - /* - for (size_t i = 0; i < ; i ++) - { - wsindexsetflag.push_back(false); - storable.push_back(false); - storL2s.push_back(0.0); - stor2Thetas.push_back(0.0); - storPhis.push_back(0.0); - //storDetids.push_back(0); - } - */ - // Map the properties from spectrum ID to workspace index for (size_t i = 0; i < specids.size(); i ++) { @@ -302,10 +308,14 @@ namespace Algorithms } // Reset detector information of each spectrum +#if 0 + NOT USED for (size_t i = 0; i < workspace->getNumberHistograms(); i ++) { if (!wsindexsetflag[i]) { + throw std::runtime_error("Ever reached?"); + // Won't be set API::ISpectrum *spectrum = workspace->getSpectrum(i); std::set<detid_t> detids = spectrum->getDetectorIDs(); @@ -347,16 +357,25 @@ namespace Algorithms } } } +#endif // Generate a new instrument // Name of the new instrument - std::string name = originstrument->getName(); - if (newinstrument) + std::string name = std::string(getProperty("InstrumentName")); + if (name.empty()) { - name = std::string(getProperty("InstrumentName")); + // Use the original L1 + if (!originstrument) + { + std::string errmsg("It is not supported that InstrumentName is not given, ", + "while there is no instrument associated to input workspace."); + g_log.error(errmsg); + throw std::runtime_error(errmsg); + } + name = originstrument->getName(); } - // b) Create a new instrument from scratch any way. + // Create a new instrument from scratch any way. Geometry::Instrument_sptr instrument(new Geometry::Instrument(name)); if (!bool(instrument)) { diff --git a/Code/Mantid/Framework/Algorithms/test/EditInstrumentGeometryTest.h b/Code/Mantid/Framework/Algorithms/test/EditInstrumentGeometryTest.h index 835106e389d6b32aa250098eba8b91c60052dc57..201347fa8bca2b27d0c3b17d4f3a848cf278046a 100644 --- a/Code/Mantid/Framework/Algorithms/test/EditInstrumentGeometryTest.h +++ b/Code/Mantid/Framework/Algorithms/test/EditInstrumentGeometryTest.h @@ -51,7 +51,6 @@ public: TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("L2","3.45") ); TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("Polar","90.09") ); TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("Azimuthal","1.84") ); - TS_ASSERT_THROWS_NOTHING( editdetector.setProperty("NewInstrument", false) ); // 4. Run TS_ASSERT_THROWS_NOTHING( editdetector.execute() ); @@ -81,29 +80,6 @@ public: } - /** Check detector parameter - */ - void checkDetectorParameters(API::MatrixWorkspace_sptr workspace, size_t wsindex, double realr, double realtth, double realphi) - { - API::ISpectrum* spectrum1 = workspace->getSpectrum(wsindex); - Geometry::Instrument_const_sptr instrument = workspace->getInstrument(); - - std::set<detid_t> detids = spectrum1->getDetectorIDs(); - TS_ASSERT_EQUALS(detids.size(), 1); - detid_t detid = 0; - std::set<detid_t>::iterator it; - for (it = detids.begin(); it != detids.end(); ++it){ - detid = *it; - } - Geometry::IDetector_const_sptr detector = instrument->getDetector(detid); - double r, tth, phi; - detector->getPos().getSpherical(r, tth, phi); - TS_ASSERT_DELTA(r, realr, 0.000001); - TS_ASSERT_DELTA(tth, realtth, 0.000001); - TS_ASSERT_DELTA(phi, realphi, 0.000001); - - } - /** Unit test to edit instrument parameters of all spectrums (>1) */ void test_MultipleWholeSpectrumEdit() @@ -120,11 +96,9 @@ public: // 3. Set Property TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("Workspace", "inputWS2") ); // TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("SpectrumIDs","3072,19456,40960,55296,74752,93184") ); - TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("SpectrumIDs","1,2,3,4,5,6") ); TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("L2","1.1,2.2,3.3,4.4,5.5,6.6") ); TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("Polar","90.1,90.2,90.3,90.4,90.5,90.6") ); TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("Azimuthal","1,2,3,4,5,6") ); - TS_ASSERT_THROWS_NOTHING( editdetector.setProperty("NewInstrument", false) ); // 4. Run TS_ASSERT_THROWS_NOTHING( editdetector.execute() ); @@ -142,32 +116,90 @@ public: } - - /** Unit test to edit instrument parameters of all spectrums (>1) + //---------------------------------------------------------------------------------------------- + /** Unit test to edit instrument parameters of all spectrums (>1) and using new detector IDs */ void test_MultiplePartialSpectrumEdit() { - // 1. Init + // Init EditInstrumentGeometry editdetector; editdetector.initialize(); - // 2. Load Data + // Load Data DataObjects::Workspace2D_sptr workspace2d = WorkspaceCreationHelper::create2DWorkspaceWithFullInstrument(6, 100, false); - API::AnalysisDataService::Instance().add("inputWS3", workspace2d); + API::AnalysisDataService::Instance().addOrReplace("inputWS3", workspace2d); - // 3.1 Set Property + // Set Property TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("Workspace", "inputWS3") ); - TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("SpectrumIDs","1,2,3") ); - TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("L2","1.1,2.2,3.3") ); - TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("Polar","90.1,90.2,90.3") ); - TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("Azimuthal","1,2,3") ); - TS_ASSERT_THROWS_NOTHING( editdetector.setProperty("NewInstrument", false) ); + TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("SpectrumIDs","1,2,3,4,5,6") ); + TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("L2","1.1,2.2,3.3, 4.4, 5.5, 6.6") ); + TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("Polar","90.1,90.2,90.3, 90.4, 90.5, 90.6") ); + TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("Azimuthal","1,2,3,4,5,6") ); + TS_ASSERT_THROWS_NOTHING( editdetector.setPropertyValue("DetectorIDs", "200, 201, 300, 301, 400, 401") ); - // 3.2 Run + // Run editdetector.execute(); - TS_ASSERT(!editdetector.isExecuted() ); + TS_ASSERT( editdetector.isExecuted() ); + + // Check + checkDetectorID(workspace2d, 0, 200); + checkDetectorID(workspace2d, 1, 201); + checkDetectorID(workspace2d, 2, 300); + checkDetectorID(workspace2d, 3, 301); + checkDetectorID(workspace2d, 4, 400); + checkDetectorID(workspace2d, 5, 401); + + // Clean + AnalysisDataService::Instance().remove("inputWS3"); + + } + + //---------------------------------------------------------------------------------------------- + /** Check detector parameter + */ + void checkDetectorParameters(API::MatrixWorkspace_sptr workspace, size_t wsindex, double realr, + double realtth, double realphi) + { + API::ISpectrum* spectrum1 = workspace->getSpectrum(wsindex); + Geometry::Instrument_const_sptr instrument = workspace->getInstrument(); + + std::set<detid_t> detids = spectrum1->getDetectorIDs(); + TS_ASSERT_EQUALS(detids.size(), 1); + detid_t detid = 0; + std::set<detid_t>::iterator it; + for (it = detids.begin(); it != detids.end(); ++it){ + detid = *it; + } + Geometry::IDetector_const_sptr detector = instrument->getDetector(detid); + double r, tth, phi; + detector->getPos().getSpherical(r, tth, phi); + TS_ASSERT_DELTA(r, realr, 0.000001); + TS_ASSERT_DELTA(tth, realtth, 0.000001); + TS_ASSERT_DELTA(phi, realphi, 0.000001); + + } + + //---------------------------------------------------------------------------------------------- + /** Check detector parameter + */ + void checkDetectorID(API::MatrixWorkspace_sptr workspace, size_t wsindex, detid_t detid) + { + API::ISpectrum* spectrum1 = workspace->getSpectrum(wsindex); + Geometry::Instrument_const_sptr instrument = workspace->getInstrument(); + + std::set<detid_t> detids = spectrum1->getDetectorIDs(); + TS_ASSERT_EQUALS(detids.size(), 1); + detid_t thisdetid = 0; + std::set<detid_t>::iterator it; + for (it = detids.begin(); it != detids.end(); ++it) + { + thisdetid = *it; + } + + TS_ASSERT_EQUALS(detid, thisdetid); + return; } };