Split yspace inputs in backward and forward

src/mvesuvio/analysis_fitting.py

@@ -15,39 +15,39 @@
 repoPath = Path(__file__).absolute().parent  # Path to the repository
 
 
-def fitInYSpaceProcedure(yFitIC, IC, wsTOF):
+def fitInYSpaceProcedure(IC, wsTOF):
 
     try:
         wsTOF = mtd[wsTOF]
     except KeyError: 
         print(f"Workspace to fit {wsTOF} not found.")
         return
 
-    wsResSum, wsRes = calculateMantidResolutionFirstMass(IC, yFitIC, wsTOF)
+    wsResSum, wsRes = calculateMantidResolutionFirstMass(IC, wsTOF)
 
     wsTOFMass0 = subtract_profiles_except_lightest(IC, wsTOF)
 
-    if yFitIC.subtract_calculated_fse_from_data:
+    if IC.subtract_calculated_fse_from_data:
         wsFSEMass0 = find_ws_name_fse_first_mass(IC)
         wsTOFMass0 = Minus(wsTOFMass0, wsFSEMass0, OutputWorkspace=wsTOFMass0.name() + "_fse")
 
-    wsJoY, wsJoYAvg = ySpaceReduction(wsTOFMass0, IC.masses[0], yFitIC, IC)
+    wsJoY, wsJoYAvg = ySpaceReduction(wsTOFMass0, IC)
 
-    if yFitIC.do_symmetrisation:
+    if IC.do_symmetrisation:
         wsJoYAvg = symmetrizeWs(wsJoYAvg)
 
-    fitProfileMinuit(yFitIC, wsJoYAvg, wsResSum)
-    fitProfileMantidFit(yFitIC, wsJoYAvg, wsResSum)
+    fitProfileMinuit(IC, wsJoYAvg, wsResSum)
+    fitProfileMantidFit(IC, wsJoYAvg, wsResSum)
 
     printYSpaceFitResults(wsJoYAvg.name())
 
-    yfitResults = ResultsYFitObject(IC, yFitIC, wsTOF.name(), wsJoYAvg.name())
+    yfitResults = ResultsYFitObject(IC, wsTOF.name(), wsJoYAvg.name())
     yfitResults.save()
 
-    if yFitIC.do_global_fit:
-        runGlobalFit(wsJoY, wsRes, IC, yFitIC)
+    if IC.do_global_fit:
+        runGlobalFit(wsJoY, wsRes, IC)
 
-    save_workspaces(yFitIC)
+    save_workspaces(IC)
     return yfitResults
 
 
@@ -67,7 +67,7 @@ def find_ws_name_fse_first_mass(ic):
     return ws_names_fse[np.argmin(ws_masses)]
 
 
-def calculateMantidResolutionFirstMass(IC, yFitIC, ws):
+def calculateMantidResolutionFirstMass(IC, ws):
     mass = IC.masses[0]
 
     resName = ws.name() + "_Resolution"
@@ -77,7 +77,7 @@ def calculateMantidResolutionFirstMass(IC, yFitIC, ws):
         )
         Rebin(
             InputWorkspace="tmp",
-            Params=yFitIC.range_for_rebinning_in_y_space,
+            Params=IC.range_for_rebinning_in_y_space,
             OutputWorkspace="tmp",
         )
 
@@ -124,19 +124,20 @@ def switchFirstTwoAxis(A):
     return np.stack(np.split(A, len(A), axis=0), axis=2)[0]
 
 
-def ySpaceReduction(wsTOF, mass0, yFitIC, ic):
+def ySpaceReduction(wsTOF, ic):
     """Seperate procedures depending on masking specified."""
 
+    mass0 = ic.masses[0]
     ws_name_lightest_mass = ic.name + '_' + str(ic.number_of_iterations_for_corrections) + '_' + str(min(ic.masses)) + '_ncp'
     ncp = mtd[ws_name_lightest_mass].extractY()
 
-    rebinPars = yFitIC.range_for_rebinning_in_y_space
+    rebinPars = ic.range_for_rebinning_in_y_space
 
     if np.any(np.all(wsTOF.extractY() == 0, axis=0)):  # Masked columns present
-        if yFitIC.mask_zeros_with == "nan":
+        if ic.mask_zeros_with == "nan":
             # Build special workspace to store accumulated points
             wsJoY = convertToYSpace(wsTOF, mass0)
-            xp = buildXRangeFromRebinPars(yFitIC)
+            xp = buildXRangeFromRebinPars(ic)
             wsJoYB = dataXBining(
                 wsJoY, xp
             )  # Unusual ws with several dataY points per each dataX point
@@ -1319,7 +1320,7 @@ def printYSpaceFitResults(wsJoYName):
 
 
 class ResultsYFitObject:
-    def __init__(self, ic, yFitIC, wsFinalName, wsYSpaceAvgName):
+    def __init__(self, ic, wsFinalName, wsYSpaceAvgName):
         # Extract most relevant information from ws
         wsFinal = mtd[wsFinalName]
         wsResSum = mtd[wsFinalName + "_Resolution_Sum"]
@@ -1339,19 +1340,19 @@ def __init__(self, ic, yFitIC, wsFinalName, wsYSpaceAvgName):
         poptList = []
         perrList = []
         try:
-            wsFitMinuit = mtd[wsJoYAvg.name() + "_minuit_" + yFitIC.fitting_model + "_Parameters"]
+            wsFitMinuit = mtd[wsJoYAvg.name() + "_minuit_" + ic.fitting_model + "_Parameters"]
             poptList.append(wsFitMinuit.column("Value"))
             perrList.append(wsFitMinuit.column("Error"))
         except:
             pass
         try:
-            wsFitLM = mtd[wsJoYAvg.name() + "_lm_" + yFitIC.fitting_model + "_Parameters"]
+            wsFitLM = mtd[wsJoYAvg.name() + "_lm_" + ic.fitting_model + "_Parameters"]
             poptList.append(wsFitLM.column("Value"))
             perrList.append(wsFitLM.column("Error"))
         except:
             pass
         try:
-            wsFitSimplex = mtd[wsJoYAvg.name() + "_simplex_" + yFitIC.fitting_model + "_Parameters"]
+            wsFitSimplex = mtd[wsJoYAvg.name() + "_simplex_" + ic.fitting_model + "_Parameters"]
             poptList.append(wsFitSimplex.column("Value"))
             perrList.append(wsFitSimplex.column("Error"))
         except:
@@ -1371,7 +1372,7 @@ def __init__(self, ic, yFitIC, wsFinalName, wsYSpaceAvgName):
         self.perr = perr
 
         self.savePath = ic.ySpaceFitSavePath
-        self.fitting_model= yFitIC.fitting_model
+        self.fitting_model= ic.fitting_model
 
     def save(self):
         np.savez(
@@ -1387,23 +1388,23 @@ def save(self):
         )
 
 
-def runGlobalFit(wsYSpace, wsRes, IC, yFitIC):
+def runGlobalFit(wsYSpace, wsRes, IC):
     print("\nRunning GLobal Fit ...\n")
 
     dataX, dataY, dataE, dataRes, instrPars = extractData(wsYSpace, wsRes, IC)
     dataX, dataY, dataE, dataRes, instrPars = takeOutMaskedSpectra(
         dataX, dataY, dataE, dataRes, instrPars
     )
 
-    idxList = groupDetectors(instrPars, yFitIC)
+    idxList = groupDetectors(instrPars, IC)
     dataX, dataY, dataE, dataRes = avgWeightDetGroups(
-        dataX, dataY, dataE, dataRes, idxList, yFitIC
+        dataX, dataY, dataE, dataRes, idxList, IC
     )
 
-    if yFitIC.do_symmetrisation:
+    if IC.do_symmetrisation:
         dataY, dataE = weightedSymArr(dataY, dataE)
 
-    model, defaultPars, sharedPars = selectModelAndPars(yFitIC.fitting_model)
+    model, defaultPars, sharedPars = selectModelAndPars(IC.fitting_model)
 
     totCost = 0
     for i, (x, y, yerr, res) in enumerate(zip(dataX, dataY, dataE, dataRes)):
@@ -1424,12 +1425,12 @@ def runGlobalFit(wsYSpace, wsRes, IC, yFitIC):
     for i in range(len(dataY)):  # Set limits for unshared parameters
         m.limits["A" + str(i)] = (0, np.inf)
 
-    if yFitIC.fitting_model== "doublewell":
+    if IC.fitting_model== "doublewell":
         m.limits["d"] = (0, np.inf)  # Shared parameters
         m.limits["R"] = (0, np.inf)
 
     t0 = time.time()
-    if yFitIC.fitting_model== "gauss":
+    if IC.fitting_model== "gauss":
         m.simplex()
         m.migrad()
 
@@ -1477,8 +1478,8 @@ def constr(*pars):
 
     create_table_for_global_fit_parameters(wsYSpace.name(), m, chi2)
 
-    if yFitIC.show_plots:
-        plotGlobalFit(dataX, dataY, dataE, m, totCost, wsYSpace.name(), yFitIC)
+    if IC.show_plots:
+        plotGlobalFit(dataX, dataY, dataE, m, totCost, wsYSpace.name(), IC)
 
     # Pass into array to store values in variable
     return np.array(m.values), np.array(m.errors)  

src/mvesuvio/config/analysis_inputs.py

@@ -11,8 +11,8 @@ class SampleParameters:
 
 @dataclass
 class BackwardAnalysisInputs(SampleParameters):
-    run_this_scattering_type = False
-    fit_in_y_space = False
+    run_this_scattering_type = True
+    fit_in_y_space = True
 
     runs = "43066-43076"
     empty_runs = "41876-41923" 
@@ -44,6 +44,28 @@ class BackwardAnalysisInputs(SampleParameters):
     multiple_scattering_number_of_events = 1.0e5
     do_gamma_correction = False
 
+    show_plots = True
+    do_symmetrisation = False
+    subtract_calculated_fse_from_data = True
+    range_for_rebinning_in_y_space = "-25, 0.5, 25"  # Needs to be symetric
+    # Fitting model options
+    # 'gauss': Single Gaussian
+    # 'gcc4': Gram-Charlier with C4 parameter
+    # 'gcc6': Gram-Charlier with C6 parameter
+    # 'doublewell': Double Well function 
+    # 'gauss3D': 3-Dimensional Gaussian 
+    fitting_model = "gauss"
+    run_minos = True
+    do_global_fit = True   # Performs global fit with Minuit by default
+    # Number of groups of detectors to perform global (simultaneous) fit on
+    # Either an integer less than the number of detectors 
+    # or option 'all', which does not form groups and fits all spectra simultaneously and individualy 
+    number_of_global_fit_groups = 4
+    # Type of masking 
+    # 'nan': Zeros in workspace being fit are ignored
+    # 'ncp': Zeros in workspace being fit are replaced by the fitted neutron compton profile 
+    mask_zeros_with = "nan"   
+
 
 @dataclass
 class ForwardAnalysisInputs(SampleParameters):
@@ -80,9 +102,6 @@ class ForwardAnalysisInputs(SampleParameters):
     multiple_scattering_number_of_events = 1.0e5
     do_gamma_correction = True
 
-
-@dataclass
-class YSpaceFitInputs:
     show_plots = True
     do_symmetrisation = False
     subtract_calculated_fse_from_data = True

src/mvesuvio/main/run_routine.py

@@ -32,7 +32,6 @@ def setup(self):
 
         self.bckwd_ai = ai.BackwardAnalysisInputs
         self.fwd_ai = ai.ForwardAnalysisInputs
-        self.yFitIC = ai.YSpaceFitInputs
 
         # Names of workspaces to check if they exist to skip analysis
         self.ws_to_fit_y_space = []
@@ -63,7 +62,8 @@ def setup(self):
         # TODO: sort out yfit inputs
         figSavePath = self.experiment_path / "figures"
         figSavePath.mkdir(exist_ok=True)
-        self.yFitIC.figSavePath = figSavePath
+        self.fwd_ai.figSavePath = figSavePath
+        self.bckwd_ai.figSavePath = figSavePath
 
 
     def import_from_inputs(self):
@@ -95,7 +95,7 @@ def run(self):
 
     def runAnalysisFitting(self):
         for wsName, i_cls in zip(self.ws_to_fit_y_space, self.classes_to_fit_y_space):
-            self.fitting_result = fitInYSpaceProcedure(self.yFitIC, i_cls, wsName)
+            self.fitting_result = fitInYSpaceProcedure(i_cls, wsName)
         return
 
 

tests/data/analysis/inputs/analysis_test.py

@@ -42,6 +42,17 @@ class BackwardAnalysisInputs(SampleParameters):
     do_gamma_correction = False
 
 
+    show_plots = False
+    do_symmetrisation = True
+    subtract_calculated_fse_from_data = False
+    range_for_rebinning_in_y_space = "-20, 0.5, 20"  # Needs to be symetric
+    fitting_model = "SINGLE_GAUSSIAN"
+    run_minos = True
+    do_global_fit = None
+    number_of_global_fit_groups = 4
+    mask_zeros_with = None
+
+
 class ForwardAnalysisInputs(SampleParameters):
     run_this_scattering_type = True
     fit_in_y_space = False
@@ -82,7 +93,6 @@ class ForwardAnalysisInputs(SampleParameters):
     do_gamma_correction= True
 
 
-class YSpaceFitInputs:
     show_plots = False
     do_symmetrisation = True
     subtract_calculated_fse_from_data = False

tests/data/analysis/inputs/yspace_gauss_test.py

@@ -1,13 +1,13 @@
 from tests.data.analysis.inputs.analysis_test import (
     BackwardAnalysisInputs,
     ForwardAnalysisInputs,
-    YSpaceFitInputs,
 )
 ForwardAnalysisInputs.fit_in_y_space = True 
 BackwardAnalysisInputs.fit_in_y_space = False
 ForwardAnalysisInputs.run_this_scattering_type = True 
 BackwardAnalysisInputs.run_this_scattering_type = False
-YSpaceFitInputs.fitting_model = "gauss"
+ForwardAnalysisInputs.fitting_model = "gauss"
+BackwardAnalysisInputs.fitting_model = "gauss"
 
 if (__name__ == "__main__") or (__name__ == "mantidqt.widgets.codeeditor.execution"):
     import mvesuvio

tests/data/analysis/inputs/yspace_gc_test.py

@@ -1,15 +1,16 @@
 from tests.data.analysis.inputs.analysis_test import (
     BackwardAnalysisInputs,
     ForwardAnalysisInputs,
-    YSpaceFitInputs,
 )
 
 ForwardAnalysisInputs.fit_in_y_space = True 
 BackwardAnalysisInputs.fit_in_y_space = False
 ForwardAnalysisInputs.run_this_scattering_type = True 
 BackwardAnalysisInputs.run_this_scattering_type = False
-YSpaceFitInputs.fitting_model = "gcc4c6"
-YSpaceFitInputs.do_symmetrisation = False 
+ForwardAnalysisInputs.fitting_model = "gcc4c6"
+ForwardAnalysisInputs.do_symmetrisation = False 
+BackwardAnalysisInputs.fitting_model = "gcc4c6"
+BackwardAnalysisInputs.do_symmetrisation = False 
 
 
 if (__name__ == "__main__") or (__name__ == "mantidqt.widgets.codeeditor.execution"):