diff --git a/Framework/WorkflowAlgorithms/inc/MantidWorkflowAlgorithms/MuonCalculateAsymmetry.h b/Framework/WorkflowAlgorithms/inc/MantidWorkflowAlgorithms/MuonCalculateAsymmetry.h index 70f448d5433155eb5d1fd7ca0cb2d466a832e7bb..0d8c35783d1c54ea3aad93fc4231da183b59000b 100644 --- a/Framework/WorkflowAlgorithms/inc/MantidWorkflowAlgorithms/MuonCalculateAsymmetry.h +++ b/Framework/WorkflowAlgorithms/inc/MantidWorkflowAlgorithms/MuonCalculateAsymmetry.h @@ -49,12 +49,22 @@ private: void init(); void exec(); - /// Converts given workspace according to the OutputType - API::MatrixWorkspace_sptr convertWorkspace(API::MatrixWorkspace_sptr ws); - - /// Merges two period workspaces according to PeriodOperation specified - API::MatrixWorkspace_sptr mergePeriods(API::MatrixWorkspace_sptr ws1, - API::MatrixWorkspace_sptr ws2); + // Calculates raw counts + API::MatrixWorkspace_sptr + calculateGroupCounts(const API::MatrixWorkspace_sptr &firstPeriodWS, + const API::MatrixWorkspace_sptr &secondPeriodWS, + int groupIndex, std::string op); + // Calculates asymmetry for specified spectrum + API::MatrixWorkspace_sptr + calculateGroupAsymmetry(const API::MatrixWorkspace_sptr &firstPeriodWS, + const API::MatrixWorkspace_sptr &secondPeriodWS, + int groupIndex, std::string op); + // Calculates asymmetry for a pair of spectra + API::MatrixWorkspace_sptr + calculatePairAsymmetry(const API::MatrixWorkspace_sptr &firstPeriodWS, + const API::MatrixWorkspace_sptr &secondPeriodWS, + int firstPairIndex, int secondPairIndex, double alpha, + std::string op); }; } // namespace WorkflowAlgorithms diff --git a/Framework/WorkflowAlgorithms/src/MuonCalculateAsymmetry.cpp b/Framework/WorkflowAlgorithms/src/MuonCalculateAsymmetry.cpp index 1384d1b189e7e8c8515a2349575377c363926b1d..66e012d86d4b76abedc22d7cf9b5e7a4b5dafe3e 100644 --- a/Framework/WorkflowAlgorithms/src/MuonCalculateAsymmetry.cpp +++ b/Framework/WorkflowAlgorithms/src/MuonCalculateAsymmetry.cpp @@ -107,126 +107,272 @@ void MuonCalculateAsymmetry::exec() { MatrixWorkspace_sptr firstPeriodWS = getProperty("FirstPeriodWorkspace"); MatrixWorkspace_sptr secondPeriodWS = getProperty("SecondPeriodWorkspace"); - MatrixWorkspace_sptr firstConverted = convertWorkspace(firstPeriodWS); + // The type of calculation + const std::string type = getPropertyValue("OutputType"); - if (secondPeriodWS) { - // Two periods - MatrixWorkspace_sptr secondConverted = convertWorkspace(secondPeriodWS); + // The group index + int groupIndex = getProperty("GroupIndex"); + + // The type of period operation (+ or -) + std::string op = getProperty("PeriodOperation"); + + if (type == "GroupCounts") { - setProperty("OutputWorkspace", - mergePeriods(firstConverted, secondConverted)); + auto outWS = + calculateGroupCounts(firstPeriodWS, secondPeriodWS, groupIndex, op); + + setProperty("OutputWorkspace", outWS); + + } else if (type == "GroupAsymmetry") { + + auto outWS = + calculateGroupAsymmetry(firstPeriodWS, secondPeriodWS, groupIndex, op); + + setProperty("OutputWorkspace", outWS); + + } else if (type == "PairAsymmetry") { + + int pairFirstIndex = getProperty("PairFirstIndex"); + int pairSecondIndex = getProperty("PairSecondIndex"); + double alpha = getProperty("Alpha"); + + auto outWS = + calculatePairAsymmetry(firstPeriodWS, secondPeriodWS, pairFirstIndex, + pairSecondIndex, alpha, op); + + setProperty("OutputWorkspace", outWS); } else { - // Single period only - setProperty("OutputWorkspace", firstConverted); - } + throw std::invalid_argument("Specified OutputType is not supported"); + } } /** - * Converts given workspace according to the OutputType. - * @param ws :: Workspace to convert - * @return Converted workspace - */ -MatrixWorkspace_sptr -MuonCalculateAsymmetry::convertWorkspace(MatrixWorkspace_sptr ws) { - const std::string type = getPropertyValue("OutputType"); +* Calculates raw counts according to period operation +* @param firstPeriodWS :: [input] First period workspace +* @param secondPeriodWS :: [input] Second period workspace +* @param groupIndex :: [input] Index of the workspace to extract counts from +* @param op :: [input] Period operation (+ or -) +*/ +MatrixWorkspace_sptr MuonCalculateAsymmetry::calculateGroupCounts( + const MatrixWorkspace_sptr &firstPeriodWS, + const MatrixWorkspace_sptr &secondPeriodWS, int groupIndex, + std::string op) { + + if (secondPeriodWS) { + // Two periods supplied + + MatrixWorkspace_sptr tempWS; - if (type == "GroupCounts" || type == "GroupAsymmetry") { - int groupIndex = getProperty("GroupIndex"); + if (op == "+") { - if (groupIndex == EMPTY_INT()) - throw std::runtime_error("GroupIndex is not specified"); + IAlgorithm_sptr alg = createChildAlgorithm("Plus"); + alg->initialize(); + alg->setProperty("LHSWorkspace", firstPeriodWS); + alg->setProperty("RHSWorkspace", secondPeriodWS); + alg->execute(); + tempWS = alg->getProperty("OutputWorkspace"); + + } else if (op == "-") { + + IAlgorithm_sptr alg = createChildAlgorithm("Minus"); + alg->initialize(); + alg->setProperty("LHSWorkspace", firstPeriodWS); + alg->setProperty("RHSWorkspace", secondPeriodWS); + alg->execute(); + tempWS = alg->getProperty("OutputWorkspace"); + } - // Yank out the counts of requested group IAlgorithm_sptr alg = createChildAlgorithm("ExtractSingleSpectrum"); alg->initialize(); - alg->setProperty("InputWorkspace", ws); + alg->setProperty("InputWorkspace", tempWS); alg->setProperty("WorkspaceIndex", groupIndex); alg->execute(); + MatrixWorkspace_sptr outWS = alg->getProperty("OutputWorkspace"); + + return outWS; + } else { + // Only one period supplied + + IAlgorithm_sptr alg = createChildAlgorithm("ExtractSingleSpectrum"); + alg->initialize(); + alg->setProperty("InputWorkspace", firstPeriodWS); + alg->setProperty("WorkspaceIndex", groupIndex); + alg->execute(); MatrixWorkspace_sptr outWS = alg->getProperty("OutputWorkspace"); - if (type == "GroupAsymmetry") { - // GroupAsymmetry - counts with ExpDecay removed and normalized + return outWS; + } +} + +/** +* Calculates single-spectrum asymmetry according to period operation +* @param firstPeriodWS :: [input] First period workspace +* @param secondPeriodWS :: [input] Second period workspace +* @param groupIndex :: [input] Workspace index for which to calculate asymmetry +* @param op :: [input] Period operation (+ or -) +*/ +MatrixWorkspace_sptr MuonCalculateAsymmetry::calculateGroupAsymmetry( + const MatrixWorkspace_sptr &firstPeriodWS, + const MatrixWorkspace_sptr &secondPeriodWS, int groupIndex, + std::string op) { + + // The output workspace + MatrixWorkspace_sptr tempWS; + + if (secondPeriodWS) { + // Two period workspaces where supplied + + if (op == "+") { + + // Sum + IAlgorithm_sptr alg = createChildAlgorithm("Plus"); + alg->initialize(); + alg->setProperty("LHSWorkspace", firstPeriodWS); + alg->setProperty("RHSWorkspace", secondPeriodWS); + alg->execute(); + MatrixWorkspace_sptr sumWS = alg->getProperty("OutputWorkspace"); + + // Remove decay + IAlgorithm_sptr asym = createChildAlgorithm("RemoveExpDecay"); + asym->setProperty("InputWorkspace",sumWS); + asym->execute(); + tempWS = asym->getProperty("OutputWorkspace"); + + } else if (op == "-") { + + // Remove decay (first period ws) + IAlgorithm_sptr asym = createChildAlgorithm("RemoveExpDecay"); + asym->initialize(); + asym->setProperty("InputWorkspace",firstPeriodWS); + asym->execute(); + MatrixWorkspace_sptr asymFirstPeriod = asym->getProperty("OutputWorkspace"); + + // Remove decay (second period ws) + asym = createChildAlgorithm("RemoveExpDecay"); + asym->initialize(); + asym->setProperty("InputWorkspace", secondPeriodWS); + asym->execute(); + MatrixWorkspace_sptr asymSecondPeriod = asym->getProperty("OutputWorkspace"); + + // Now subtract + IAlgorithm_sptr alg = createChildAlgorithm("Minus"); + alg->initialize(); + alg->setProperty("LHSWorkspace", asymFirstPeriod); + alg->setProperty("RHSWorkspace", asymSecondPeriod); + alg->execute(); + tempWS = alg->getProperty("OutputWorkspace"); + + } + } else { + // Only one period was supplied IAlgorithm_sptr alg = createChildAlgorithm("RemoveExpDecay"); alg->initialize(); - alg->setProperty("InputWorkspace", outWS); + alg->setProperty("InputWorkspace", firstPeriodWS); alg->execute(); - - outWS = alg->getProperty("OutputWorkspace"); + tempWS = alg->getProperty("OutputWorkspace"); } + // Extract the requested spectrum + IAlgorithm_sptr alg = createChildAlgorithm("ExtractSingleSpectrum"); + alg->initialize(); + alg->setProperty("InputWorkspace", tempWS); + alg->setProperty("WorkspaceIndex", groupIndex); + alg->execute(); + MatrixWorkspace_sptr outWS = alg->getProperty("OutputWorkspace"); return outWS; - } else if (type == "PairAsymmetry") { - // PairAsymmetry - result of AsymmetryCalc algorithm +} - int pairFirstIndex = getProperty("PairFirstIndex"); - int pairSecondIndex = getProperty("PairSecondIndex"); +/** +* Calculates pair asymmetry according to period operation +* @param firstPeriodWS :: [input] First period workspace +* @param secondPeriodWS :: [input] Second period workspace +* @param groupIndex :: [input] Workspace index for which to calculate asymmetry +* @param op :: [input] Period operation (+ or -) +*/ +MatrixWorkspace_sptr MuonCalculateAsymmetry::calculatePairAsymmetry( + const MatrixWorkspace_sptr &firstPeriodWS, + const MatrixWorkspace_sptr &secondPeriodWS, int firstPairIndex, + int secondPairIndex, double alpha, std::string op) { + + // Pair indices as vectors + std::vector<int> fwd(1, firstPairIndex + 1); + std::vector<int> bwd(1, secondPairIndex + 1); - if (pairFirstIndex == EMPTY_INT() || pairSecondIndex == EMPTY_INT()) - throw std::invalid_argument("Both pair indices should be specified"); + if (secondPeriodWS) { - double alpha = getProperty("Alpha"); + MatrixWorkspace_sptr outWS; + + if (op == "+") { + + // Sum + IAlgorithm_sptr alg = createChildAlgorithm("Plus"); + alg->initialize(); + alg->setProperty("LHSWorkspace", firstPeriodWS); + alg->setProperty("RHSWorkspace", secondPeriodWS); + alg->execute(); + MatrixWorkspace_sptr sumWS = alg->getProperty("OutputWorkspace"); + + // Asymmetry calculation + alg = createChildAlgorithm("AsymmetryCalc"); + alg->setProperty("InputWorkspace", sumWS); + alg->setProperty("ForwardSpectra", fwd); + alg->setProperty("BackwardSpectra", bwd); + alg->setProperty("Alpha", alpha); + alg->execute(); + outWS = alg->getProperty("OutputWorkspace"); + + } else if (op == "-") { + + std::vector<int> fwd(1, firstPairIndex + 1); + std::vector<int> bwd(1, secondPairIndex + 1); + + // First period asymmetry + IAlgorithm_sptr alg = createChildAlgorithm("AsymmetryCalc"); + alg->setProperty("InputWorkspace", firstPeriodWS); + alg->setProperty("ForwardSpectra", fwd); + alg->setProperty("BackwardSpectra", bwd); + alg->setProperty("Alpha", alpha); + alg->execute(); + MatrixWorkspace_sptr asymFirstPeriod = + alg->getProperty("OutputWorkspace"); + + // Second period asymmetry + alg = createChildAlgorithm("AsymmetryCalc"); + alg->setProperty("InputWorkspace", secondPeriodWS); + alg->setProperty("ForwardSpectra", fwd); + alg->setProperty("BackwardSpectra", bwd); + alg->setProperty("Alpha", alpha); + alg->execute(); + MatrixWorkspace_sptr asymSecondPeriod = + alg->getProperty("OutputWorkspace"); + + // Now subtract + alg = createChildAlgorithm("Minus"); + alg->initialize(); + alg->setProperty("LHSWorkspace", asymFirstPeriod); + alg->setProperty("RHSWorkspace", asymSecondPeriod); + alg->execute(); + outWS = alg->getProperty("OutputWorkspace"); + } - // We get pair groups as their workspace indices, but AsymmetryCalc wants - // spectra numbers, - // so need to convert - specid_t spectraNo1 = ws->getSpectrum(pairFirstIndex)->getSpectrumNo(); - specid_t spectraNo2 = ws->getSpectrum(pairSecondIndex)->getSpectrumNo(); + return outWS; - if (spectraNo1 == -1 || spectraNo2 == -1 || spectraNo1 == spectraNo2) - throw std::invalid_argument( - "Spectra numbers of the input workspace are not set properly"); + } else { IAlgorithm_sptr alg = createChildAlgorithm("AsymmetryCalc"); - alg->setProperty("InputWorkspace", ws); - // As strings, cause otherwise would need to create arrays with single - // elements - alg->setPropertyValue("ForwardSpectra", - boost::lexical_cast<std::string>(spectraNo1)); - alg->setPropertyValue("BackwardSpectra", - boost::lexical_cast<std::string>(spectraNo2)); + alg->setProperty("InputWorkspace", firstPeriodWS); + alg->setProperty("ForwardSpectra", fwd); + alg->setProperty("BackwardSpectra", bwd); alg->setProperty("Alpha", alpha); alg->execute(); - MatrixWorkspace_sptr outWS = alg->getProperty("OutputWorkspace"); return outWS; } - - throw std::invalid_argument("Specified OutputType is not supported"); -} - -/** - * Merges two period workspaces according to PeriodOperation specified. - * @param ws1 :: First period workspace - * @param ws2 :: Second period workspace - * @return Merged workspace - */ -MatrixWorkspace_sptr -MuonCalculateAsymmetry::mergePeriods(MatrixWorkspace_sptr ws1, - MatrixWorkspace_sptr ws2) { - std::string op = getProperty("PeriodOperation"); - - std::string algorithmName; - - if (op == "+") { - algorithmName = "Plus"; - } else if (op == "-") { - algorithmName = "Minus"; - } - - IAlgorithm_sptr alg = createChildAlgorithm(algorithmName); - alg->initialize(); - alg->setProperty("LHSWorkspace", ws1); - alg->setProperty("RHSWorkspace", ws2); - alg->execute(); - - MatrixWorkspace_sptr outWS = alg->getProperty("OutputWorkspace"); - - return outWS; } } // namespace WorkflowAlgorithms } // namespace Mantid diff --git a/Framework/WorkflowAlgorithms/test/MuonCalculateAsymmetryTest.h b/Framework/WorkflowAlgorithms/test/MuonCalculateAsymmetryTest.h index 9adc9e49216d95b5f11bc9b641aeddf60f39f730..57dc1e3db8fd52e17d222e4603e9b0e29891d776 100644 --- a/Framework/WorkflowAlgorithms/test/MuonCalculateAsymmetryTest.h +++ b/Framework/WorkflowAlgorithms/test/MuonCalculateAsymmetryTest.h @@ -317,6 +317,7 @@ public: alg.initialize(); alg.setProperty("FirstPeriodWorkspace", inWS); alg.setProperty("SecondPeriodWorkspace", inWSSecond); + alg.setProperty("PeriodOperation", "-"); alg.setProperty("OutputType", "PairAsymmetry"); alg.setProperty("PeriodOperation", "-"); alg.setProperty("PairFirstIndex", 2);