Put set up of fit parameters in class methods

Initialize fit parameters and bounds when setting up the class.
This makes it easier to keep track of the order in which the masses are
getting into the parameters array and makes it easier to change it in the
future.

Made self._masses method because it is used several times across the routine.

src/mvesuvio/analysis_reduction.py

@@ -127,6 +127,26 @@ def _setup(self):
 
         self._profiles_table = self.getProperty("InputProfiles").value
 
+        # Need to transform profiles table into parameter array for minimize
+        self._initial_fit_parameters = []
+        for intensity, width, center in zip(
+            self._profiles_table.column("intensity"),
+            self._profiles_table.column("width"),
+            self._profiles_table.column("center")
+        ):
+            self._initial_fit_parameters.extend([intensity, width, center])
+
+        self._initial_fit_bounds = []
+        for intensity_bounds, width_bounds, center_bounds in zip(
+            self._profiles_table.column("intensity_bounds"),
+            self._profiles_table.column("width_bounds"),
+            self._profiles_table.column("center_bounds")
+        ):
+            self._initial_fit_bounds.extend([eval(intensity_bounds), eval(width_bounds), eval(center_bounds)])
+
+        # Masses need to be defined in the same order
+        self._masses = np.array(self._profiles_table.column("mass"))
+
         # Variables changing during fit
         self._workspace_for_corrections = self.getProperty("InputWorkspace").value 
         self._workspace_being_fit = self.getProperty("InputWorkspace").value
@@ -142,16 +162,15 @@ def _setup(self):
 
     def calculate_h_ratio(self):
 
-        if not np.isclose(min(self._profiles_table.column("mass")), 1, atol=0.1):    # Hydrogen not present
+        if not np.isclose(min(self._masses), 1, atol=0.1):    # Hydrogen not present
             return None
 
         # Hydrogen present 
         # intensities = np.array([p.mean_intensity for p in self._profiles.values()])
         intensities = self.mean_intensity_ratios
         # masses = np.array([p.mass for p in self._profiles.values()])
-        masses = self._profiles_table.column("mass")
 
-        sorted_intensities = intensities[np.argsort(masses)]
+        sorted_intensities = intensities[np.argsort(self._masses)]
         return sorted_intensities[0] / sorted_intensities[1] 
 
 
@@ -367,7 +386,7 @@ def convertDataXToYSpacesForEachMass(self, dataX, delta_Q, delta_E):
         delta_E = delta_E[np.newaxis, :, :]
 
         mN, Ef, en_to_vel, vf, hbar = loadConstants()
-        masses = np.array(self._profiles_table.column("mass")).reshape(-1, 1, 1)
+        masses = self._masses.reshape(-1, 1, 1)
 
         energyRecoil = np.square(hbar * delta_Q) / 2.0 / masses
         ySpacesForEachMass = (
@@ -543,28 +562,11 @@ def _fit_neutron_compton_profiles_to_row(self):
             self._table_fit_results.addRow(np.zeros(3*self._profiles_table.rowCount()+3))
             return
 
-        # Need to transform profiles into parameter array for minimize
-        initial_parameters = []
-        for intensity, width, center in zip(
-            self._profiles_table.column("intensity"),
-            self._profiles_table.column("width"),
-            self._profiles_table.column("center")
-        ):
-            initial_parameters.extend([intensity, width, center])
-
-        bounds = []
-        for intensity_bounds, width_bounds, center_bounds in zip(
-            self._profiles_table.column("intensity_bounds"),
-            self._profiles_table.column("width_bounds"),
-            self._profiles_table.column("center_bounds")
-        ):
-            bounds.extend([eval(intensity_bounds), eval(width_bounds), eval(center_bounds)])
-
         result = optimize.minimize(
             self.errorFunction,
-            initial_parameters,
+            self._initial_fit_parameters,
             method="SLSQP",
-            bounds=bounds,
+            bounds=self._initial_fit_bounds,
             constraints=self._constraints,
         )
         fitPars = result["x"]
@@ -618,7 +620,7 @@ def _neutron_compton_profiles(self, pars):
         intensities = pars[::3].reshape(-1, 1)
         widths = pars[1::3].reshape(-1, 1)
         centers = pars[2::3].reshape(-1, 1)
-        masses = np.array(self._profiles_table.column("mass")).reshape(-1, 1)
+        masses = self._masses.reshape(-1, 1)
 
         v0, E0, deltaE, deltaQ = self._kinematic_arrays[self._row_being_fit]
 
@@ -669,7 +671,7 @@ def kinematicsAtYCenters(self, centers):
 
 
     def calcGaussianResolution(self, centers):
-        masses = np.array(self._profiles_table.column("mass")).reshape(-1, 1)
+        masses = self._masses.reshape(-1, 1)
         v0, E0, delta_E, delta_Q = self.kinematicsAtYCenters(centers)
         det, plick, angle, T0, L0, L1 = self._instrument_params[self._row_being_fit]
         dE1, dTOF, dTheta, dL0, dL1, dE1_lorz = self._resolution_params[self._row_being_fit]
@@ -715,7 +717,7 @@ def calcGaussianResolution(self, centers):
 
 
     def calcLorentzianResolution(self, centers):
-        masses = np.array(self._profiles_table.column("mass")).reshape(-1, 1)
+        masses = self._masses.reshape(-1, 1)
         v0, E0, delta_E, delta_Q = self.kinematicsAtYCenters(centers)
         det, plick, angle, T0, L0, L1 = self._instrument_params[self._row_being_fit]
         dE1, dTOF, dTheta, dL0, dL1, dE1_lorz = self._resolution_params[self._row_being_fit]
@@ -884,7 +886,7 @@ def createSlabGeometry(self, wsNCPM):
 
 
     def calcMSCorrectionSampleProperties(self, meanWidths, meanIntensityRatios):
-        masses = np.array(self._profiles_table.column("mass"))
+        masses = self._masses
 
         # If Backscattering mode and H is present in the sample, add H to MS properties
         if self._mode_running == "BACKWARD":
@@ -975,9 +977,8 @@ def create_gamma_workspaces(self):
 
 
     def calcGammaCorrectionProfiles(self, meanWidths, meanIntensityRatios):
-        masses = np.array(self._profiles_table.column("mass"))
         profiles = ""
-        for mass, width, intensity in zip(masses, meanWidths, meanIntensityRatios):
+        for mass, width, intensity in zip(self._masses, meanWidths, meanIntensityRatios):
             profiles += (
                 "name=GaussianComptonProfile,Mass="
                 + str(mass)

src/mvesuvio/analysis_routines.py

@@ -83,11 +83,10 @@ def create_profiles_table(name, profiles: list[NeutronComptonProfile]):
     table.addColumn(type="float", name="center")
     table.addColumn(type="str", name="center_bounds")
     for p in profiles:
-        table.addRow([str(getattr(p, attr)) if "bounds" in attr else getattr(p, attr) for attr in table.getColumnNames()])
-
-    for p in profiles:
-        print(str(getattr(p, "intensity_bounds")))
-        print(str(getattr(p, "width_bounds")))
+        table.addRow([str(getattr(p, attr)) 
+            if "bounds" in attr 
+            else getattr(p, attr) 
+            for attr in table.getColumnNames()])
     return table