Re #13952 Fix period addition, needed to refactor the code

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

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

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);