diff --git a/Code/Mantid/Framework/Crystal/CMakeLists.txt b/Code/Mantid/Framework/Crystal/CMakeLists.txt index 3792cfb33636..d6abd0fdfd3a 100644 --- a/Code/Mantid/Framework/Crystal/CMakeLists.txt +++ b/Code/Mantid/Framework/Crystal/CMakeLists.txt @@ -21,6 +21,7 @@ set ( SRC_FILES src/IndexPeaks.cpp src/IndexSXPeaks.cpp src/IntegratePeakTimeSlices.cpp + src/IntegratePeaksUsingClusters.cpp src/LatticeErrors.cpp src/LoadHKL.cpp src/LoadIsawPeaks.cpp @@ -83,6 +84,7 @@ set ( INC_FILES inc/MantidCrystal/IndexPeaks.h inc/MantidCrystal/IndexSXPeaks.h inc/MantidCrystal/IntegratePeakTimeSlices.h + inc/MantidCrystal/IntegratePeaksUsingClusters.h inc/MantidCrystal/LatticeErrors.h inc/MantidCrystal/LoadHKL.h inc/MantidCrystal/LoadIsawPeaks.h @@ -142,9 +144,10 @@ set ( TEST_FILES IndexPeaksTest.h IndexSXPeaksTest.h IntegratePeakTimeSlicesTest.h + IntegratePeaksUsingClustersTest.h LoadHKLTest.h LoadIsawPeaksTest.h - LoadIsawSpectrumTest.h + LoadIsawSpectrumTest.h LoadIsawUBTest.h MaskPeaksWorkspaceTest.h NormaliseVanadiumTest.h diff --git a/Code/Mantid/Framework/Crystal/inc/MantidCrystal/IntegratePeaksUsingClusters.h b/Code/Mantid/Framework/Crystal/inc/MantidCrystal/IntegratePeaksUsingClusters.h new file mode 100644 index 000000000000..87f970edd3c9 --- /dev/null +++ b/Code/Mantid/Framework/Crystal/inc/MantidCrystal/IntegratePeaksUsingClusters.h @@ -0,0 +1,56 @@ +#ifndef MANTID_CRYSTAL_INTEGRATEPEAKSUSINGCLUSTERS_H_ +#define MANTID_CRYSTAL_INTEGRATEPEAKSUSINGCLUSTERS_H_ + +#include "MantidKernel/System.h" +#include "MantidAPI/Algorithm.h" + +namespace Mantid +{ +namespace Crystal +{ + + /** IntegratePeaksUsingClusters : Uses clustering to integrate peaks. + + Copyright © 2014 ISIS Rutherford Appleton Laboratory & NScD Oak Ridge National Laboratory + + This file is part of Mantid. + + Mantid is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 3 of the License, or + (at your option) any later version. + + Mantid is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program. If not, see . + + File change history is stored at: + Code Documentation is available at: + */ + class DLLExport IntegratePeaksUsingClusters : public API::Algorithm + { + public: + IntegratePeaksUsingClusters(); + virtual ~IntegratePeaksUsingClusters(); + + virtual const std::string name() const; + virtual int version() const; + virtual const std::string category() const; + + private: + virtual void initDocs(); + void init(); + void exec(); + + + }; + + +} // namespace Crystal +} // namespace Mantid + +#endif /* MANTID_CRYSTAL_INTEGRATEPEAKSUSINGCLUSTERS_H_ */ \ No newline at end of file diff --git a/Code/Mantid/Framework/Crystal/src/DisjointElement.cpp b/Code/Mantid/Framework/Crystal/src/DisjointElement.cpp index 7bbc6cca6d6d..3e4c03848f2f 100644 --- a/Code/Mantid/Framework/Crystal/src/DisjointElement.cpp +++ b/Code/Mantid/Framework/Crystal/src/DisjointElement.cpp @@ -37,11 +37,7 @@ namespace Mantid DisjointElement::DisjointElement(const DisjointElement& other) : m_parent(other.m_parent), m_rank(other.m_rank), m_id(other.m_id) { - if (m_rank > 0) - { - throw std::logic_error( - "This is a parent node. Children cannot be copied, leading to possible inconsistencies."); - } + // Don't point to copy object as parent if copy object is it's own parent. if (other.m_parent == &other) { @@ -59,11 +55,6 @@ namespace Mantid if (this != &other) { - if (other.m_rank > 0) - { - throw std::logic_error( - "This is a parent node. Children cannot be copied, leading to possible inconsistencies."); - } m_parent = other.m_parent; m_rank = other.m_rank; @@ -183,10 +174,7 @@ namespace Mantid void DisjointElement::unionWith(DisjointElement* other) { - if (this->getId() == other->getId()) - { - throw std::logic_error("Trying to union two elements with the same Id"); - } + if (other->getRoot() != this->getRoot()) // Check sets do not already have the same root before continuing { diff --git a/Code/Mantid/Framework/Crystal/src/IntegratePeaksUsingClusters.cpp b/Code/Mantid/Framework/Crystal/src/IntegratePeaksUsingClusters.cpp new file mode 100644 index 000000000000..49268f60df4e --- /dev/null +++ b/Code/Mantid/Framework/Crystal/src/IntegratePeaksUsingClusters.cpp @@ -0,0 +1,180 @@ +/*WIKI* +Handles integration of arbitary single crystal peaks shapes. + +Uses connected component analysis to integrate peaks in an PeaksWorkspace over an MDHistoWorkspace of data. +*WIKI*/ + +#include "MantidCrystal/IntegratePeaksUsingClusters.h" +#include "MantidCrystal/ConnectedComponentLabeling.h" +#include "MantidCrystal/HardThresholdBackground.h" +#include "MantidAPI/IMDHistoWorkspace.h" +#include "MantidAPI/WorkspaceProperty.h" +#include "MantidAPI/IMDIterator.h" +#include "MantidKernel/ListValidator.h" +#include "MantidKernel/Utils.h" +#include "MantidDataObjects/PeaksWorkspace.h" +#include + +using namespace Mantid::API; +using namespace Mantid::Kernel; +using namespace Mantid::DataObjects; + +namespace +{ + class PeakBackground : public Mantid::Crystal::HardThresholdBackground + { + private: + + IPeaksWorkspace_const_sptr m_peaksWS; + const double m_radiusEstimate; + const SpecialCoordinateSystem m_mdCoordinates; + + public: + PeakBackground(IPeaksWorkspace_const_sptr peaksWS, const double& radiusEstimate, const double& thresholdSignal, const Mantid::API::MDNormalization normalisation, const SpecialCoordinateSystem coordinates) + : Mantid::Crystal::HardThresholdBackground(thresholdSignal, normalisation), m_peaksWS(peaksWS), m_radiusEstimate(radiusEstimate), m_mdCoordinates(coordinates) + { + } + + virtual bool isBackground(Mantid::API::IMDIterator* iterator) const + { + if(!HardThresholdBackground::isBackground(iterator) ) + { + for(size_t i = 0; i < m_peaksWS->getNumberPeaks(); ++i) + { + V3D coords; + if(m_mdCoordinates==QLab) + { + coords= m_peaksWS->getPeak(i).getQLabFrame(); + } + else if(m_mdCoordinates==QSample) + { + coords= m_peaksWS->getPeak(i).getQSampleFrame(); + } + else if(m_mdCoordinates==Mantid::API::HKL) + { + coords= m_peaksWS->getPeak(i).getHKL(); + } + const VMD& center = iterator->getCenter(); + V3D temp(center[0], center[1], center[2]); + if(coords.distance(temp) < m_radiusEstimate) + { + return false; + } + } + + } + return true; + } + + void configureIterator(Mantid::API::IMDIterator* const iterator) const + { + } + + virtual ~PeakBackground() + { + } + }; +} + +namespace Mantid +{ + namespace Crystal + { + + // Register the algorithm into the AlgorithmFactory + DECLARE_ALGORITHM(IntegratePeaksUsingClusters) + + + + //---------------------------------------------------------------------------------------------- + /** Constructor + */ + IntegratePeaksUsingClusters::IntegratePeaksUsingClusters() + { + } + + //---------------------------------------------------------------------------------------------- + /** Destructor + */ + IntegratePeaksUsingClusters::~IntegratePeaksUsingClusters() + { + } + + + //---------------------------------------------------------------------------------------------- + /// Algorithm's name for identification. @see Algorithm::name + const std::string IntegratePeaksUsingClusters::name() const { return "IntegratePeaksUsingClusters";}; + + /// Algorithm's version for identification. @see Algorithm::version + int IntegratePeaksUsingClusters::version() const { return 1;}; + + /// Algorithm's category for identification. @see Algorithm::category + const std::string IntegratePeaksUsingClusters::category() const { return "MDAlgorithms";} + + //---------------------------------------------------------------------------------------------- + /// Sets documentation strings for this algorithm + void IntegratePeaksUsingClusters::initDocs() + { + this->setWikiSummary("Integrate single crystal peaks using connected component analysis"); + this->setOptionalMessage(this->getWikiSummary()); + } + + //---------------------------------------------------------------------------------------------- + /** Initialize the algorithm's properties. + */ + void IntegratePeaksUsingClusters::init() + { + declareProperty(new WorkspaceProperty("InputWorkspace","",Direction::Input), "Input md workspace."); + std::vector propOptions; + propOptions.push_back("NoNormalization"); + propOptions.push_back("NumberOfEventsNormalization"); + propOptions.push_back("VolumeNormalization"); + declareProperty("Normalization", "NoNormalization",boost::make_shared(propOptions), + "Normalization to use." + ); + + declareProperty(new WorkspaceProperty("PeaksWorkspace","", Direction::Input),"A PeaksWorkspace containing the peaks to integrate."); + declareProperty(new PropertyWithValue("RadiusEstimate", 0.1, Direction::Input), "Estimate of Peak Radius. This is optional and should not affect the ouput of the integration, but it will remove noise from the resulting MD cluster image."); + declareProperty(new PropertyWithValue("Threshold", 0, Direction::Input), "Threshold"); + //declareProperty(new WorkspaceProperty("OutputWorkspace","",Direction::Output), "An output integrated peaks workspace."); + declareProperty(new WorkspaceProperty("OutputWorkspaceMD","",Direction::Output), "MDHistoWorkspace containing the clusters."); + } + + //---------------------------------------------------------------------------------------------- + /** Execute the algorithm. + */ + void IntegratePeaksUsingClusters::exec() + { + IMDHistoWorkspace_sptr mdWS = getProperty("InputWorkspace"); + PeaksWorkspace_sptr peaksWS = getProperty("PeaksWorkspace"); + + const SpecialCoordinateSystem mdCoordinates = mdWS->getSpecialCoordinateSystem(); + + + std::string strNormalization = getPropertyValue("Normalization"); + MDNormalization normalization = NoNormalization; + if(strNormalization == "NumberOfEventsNormalization") + { + normalization = NumEventsNormalization; + } + else if(strNormalization == "VolumeNormalization") + { + normalization = VolumeNormalization; + } + + const double threshold = getProperty("Threshold"); + const double radiusEstimate = getProperty("RadiusEstimate"); + PeakBackground background(peaksWS, radiusEstimate, threshold, normalization, mdCoordinates); + //HardThresholdBackground background(threshold, normalization); + + ConnectedComponentLabeling analysis; + IMDHistoWorkspace_sptr clusters = analysis.execute(mdWS, &background); + + //setProperty("OutputWorkspace", peaksWS); + setProperty("OutputWorkspaceMD", clusters); + } + + + + } // namespace Crystal +} // namespace Mantid \ No newline at end of file diff --git a/Code/Mantid/Framework/Crystal/test/IntegratePeaksUsingClustersTest.h b/Code/Mantid/Framework/Crystal/test/IntegratePeaksUsingClustersTest.h new file mode 100644 index 000000000000..7febbc1ec977 --- /dev/null +++ b/Code/Mantid/Framework/Crystal/test/IntegratePeaksUsingClustersTest.h @@ -0,0 +1,60 @@ +#ifndef MANTID_CRYSTAL_INTEGRATEPEAKSUSINGCLUSTERSTEST_H_ +#define MANTID_CRYSTAL_INTEGRATEPEAKSUSINGCLUSTERSTEST_H_ + +#include + +#include "MantidCrystal/IntegratePeaksUsingClusters.h" + +using Mantid::Crystal::IntegratePeaksUsingClusters; + +class IntegratePeaksUsingClustersTest : public CxxTest::TestSuite +{ +public: + // This pair of boilerplate methods prevent the suite being created statically + // This means the constructor isn't called when running other tests + static IntegratePeaksUsingClustersTest *createSuite() { return new IntegratePeaksUsingClustersTest(); } + static void destroySuite( IntegratePeaksUsingClustersTest *suite ) { delete suite; } + + + void test_Init() + { + IntegratePeaksUsingClusters alg; + TS_ASSERT_THROWS_NOTHING( alg.initialize() ) + TS_ASSERT( alg.isInitialized() ) + } + + void test_exec() + { + // Name of the output workspace. + std::string outWSName("IntegratePeaksUsingClustersTest_OutputWS"); + + IntegratePeaksUsingClusters alg; + TS_ASSERT_THROWS_NOTHING( alg.initialize() ) + TS_ASSERT( alg.isInitialized() ) + TS_ASSERT_THROWS_NOTHING( alg.setPropertyValue("REPLACE_PROPERTY_NAME_HERE!!!!", "value") ); + TS_ASSERT_THROWS_NOTHING( alg.setPropertyValue("OutputWorkspace", outWSName) ); + TS_ASSERT_THROWS_NOTHING( alg.execute(); ); + TS_ASSERT( alg.isExecuted() ); + + // Retrieve the workspace from data service. TODO: Change to your desired type + Workspace_sptr ws; + TS_ASSERT_THROWS_NOTHING( ws = AnalysisDataService::Instance().retrieveWS(outWSName) ); + TS_ASSERT(ws); + if (!ws) return; + + // TODO: Check the results + + // Remove workspace from the data service. + AnalysisDataService::Instance().remove(outWSName); + } + + void test_Something() + { + TSM_ASSERT( "You forgot to write a test!", 0); + } + + +}; + + +#endif /* MANTID_CRYSTAL_INTEGRATEPEAKSUSINGCLUSTERSTEST_H_ */ \ No newline at end of file