Merge branch 'corems_massspectra_importer_fix' into 'master'

Improve ReadCoreMSHDF_MassSpectrum import and export

Closes #152 and #153

See merge request mass-spectrometry/corems!122

Author: corilo

NEG_ESI_SRFA_CoreMS_prob_score.csv

@@ -52,8 +52,6 @@ def __init__(self, file_location):
 
         self.scans = list(self.h5pydata.keys())
 
-        print(self.scans)
-
     def load_raw_data(self, mass_spectrum, scan_index=0):
         """
         Load raw data into the mass spectrum object.
@@ -68,13 +66,16 @@ def load_raw_data(self, mass_spectrum, scan_index=0):
 
         scan_label = self.scans[scan_index]
 
-        mz_profile = self.h5pydata[scan_label]["raw_ms"][0]
+        # Check if the "raw_ms" group in the scan is empty
+        if self.h5pydata[scan_label]["raw_ms"].shape is not None:
+
+            mz_profile = self.h5pydata[scan_label]["raw_ms"][0]
 
-        abundance_profile = self.h5pydata[scan_label]["raw_ms"][1]
+            abundance_profile = self.h5pydata[scan_label]["raw_ms"][1]
 
-        mass_spectrum.mz_exp_profile = mz_profile
+            mass_spectrum.mz_exp_profile = mz_profile
 
-        mass_spectrum.abundance_profile = abundance_profile
+            mass_spectrum.abundance_profile = abundance_profile
 
     def get_mass_spectrum(
         self,
@@ -83,9 +84,12 @@ def get_mass_spectrum(
         auto_process=True,
         load_settings=True,
         load_raw=True,
+        load_molecular_formula=True,
     ):
         """
-        Get a mass spectrum object.
+        Instantiate a mass spectrum object from the CoreMS HDF5 file. 
+        Note that this always returns a centroid mass spectrum object; functionality for profile and
+        frequency mass spectra is not yet implemented.
 
         Parameters
         ----------
@@ -99,14 +103,23 @@ def get_mass_spectrum(
             Whether to load the settings into the mass spectrum object. Default is True.
         load_raw : bool, optional
             Whether to load the raw data into the mass spectrum object. Default is True.
+        load_molecular_formula : bool, optional
+            Whether to load the molecular formula into the mass spectrum object.
+            Default is True.
 
         Returns
         -------
         MassSpecCentroid
             The mass spectrum object.
+        
+        Raises
+        ------
+        ValueError
+            If the CoreMS file is not valid.
+            If the mass spectrum has not been processed and load_molecular_formula is True.
         """
-
-        dataframe = self.get_dataframe(scan_number, time_index=time_index)
+        scan_index = self.scans.index(str(scan_number))
+        dataframe = self.get_dataframe(scan_index, time_index=time_index)
 
         if not set(
             ["H/C", "O/C", "Heteroatom Class", "Ion Type", "Is Isotopologue"]
@@ -117,23 +130,45 @@ def get_mass_spectrum(
 
         dataframe.rename(columns=self.parameters.header_translate, inplace=True)
 
+        # Cast m/z, and 'Peak Height' to float
+        dataframe["m/z"] = dataframe["m/z"].astype(float)
+        dataframe["Peak Height"] = dataframe["Peak Height"].astype(float)
+
         polarity = dataframe["Ion Charge"].values[0]
 
-        output_parameters = self.get_output_parameters(polarity, scan_index=scan_number)
+        output_parameters = self.get_output_parameters(polarity, scan_index=scan_index)
 
         mass_spec_obj = MassSpecCentroid(
-            dataframe.to_dict(orient="list"), output_parameters
+            dataframe.to_dict(orient="list"), output_parameters, auto_process = False
         )
 
+        if auto_process:
+            # Set the settings on the mass spectrum object to relative abuncance of 0 so all peaks get added
+            mass_spec_obj.settings.noise_threshold_method = "absolute_abundance"
+            mass_spec_obj.settings.noise_threshold_absolute_abundance = 0
+            mass_spec_obj.process_mass_spec()
+
         if load_settings:
+            # Load settings into the mass spectrum object
             self.load_settings(
-                mass_spec_obj, scan_index=scan_number, time_index=time_index
+                mass_spec_obj, 
+                scan_index=scan_index, 
+                time_index=time_index
             )
 
         if load_raw:
-            self.load_raw_data(mass_spec_obj, scan_index=scan_number)
-
-        self.add_molecular_formula(mass_spec_obj, dataframe)
+            self.load_raw_data(
+                mass_spec_obj, 
+                scan_index=scan_index
+                )
+
+        if load_molecular_formula:
+            if not auto_process:
+                raise ValueError(
+                    "Can only add molecular formula if the mass spectrum has been processed"
+                )
+            else:
+                self.add_molecular_formula(mass_spec_obj, dataframe)
 
         return mass_spec_obj
 

NEG_ESI_SRFA_CoreMS_prob_score.json

@@ -306,7 +306,8 @@ def to_hdf(self):
 
                 scan_group = hdf_handle.create_group(str(self.mass_spectrum.scan_number))
 
-                if list(self.mass_spectrum.abundance_profile):
+                # If there is raw data (from profile data) save it
+                if not self.mass_spectrum.is_centroid:
 
                     mz_abun_array = empty(shape=(2, len(self.mass_spectrum.abundance_profile)))