Added cutoff_error parameter to enable Automatic Linestrength Cutoff

radis/lbl/base.py

@@ -3087,7 +3087,7 @@ def calc_linestrength_noneq(self):
 
     # %%
     def calc_emission_integral(self):
-        r"""Calculate the emission integral (in :math:`mW/sr`) of all lines in DataFrame ``df1``.
+        """Calculate the emission integral (in :math:`mW/sr`) of all lines in DataFrame ``df1``.
 
         .. math::
 
@@ -3155,8 +3155,8 @@ def calc_emission_integral(self):
         return
 
     # %%
-    def _cutoff_linestrength(self, cutoff=None):
-        """Discard linestrengths that are lower that this, to reduce
+    def _cutoff_linestrength(self, cutoff=None, cutoff_error=None):
+        r"""Discard linestrengths that are lower that this, to reduce
         calculation times. Set the number of lines cut in
         ``self._Nlines_cutoff``
 
@@ -3167,6 +3167,10 @@ def _cutoff_linestrength(self, cutoff=None):
             discard linestrengths that are lower that this, to reduce calculation
             times. If 0, no cutoff. Default 0
 
+        cutoff_error: float
+            user-inputted value of cutoff error, to keep the estimated error
+            below this value. If None, no error to consider. Default None.
+
         Notes
         -----
 
@@ -3179,11 +3183,17 @@ def _cutoff_linestrength(self, cutoff=None):
         # Update defaults
         if cutoff is not None:
             self.params.cutoff = cutoff
+        if cutoff_error is not None:
+            self.params.cutoff_error = cutoff_error
 
         # Load variables
         cutoff = self.params.cutoff
+        cutoff_error = self.params.cutoff_error
         verbose = self.verbose
         df = self.df1
+        df1 = self.df1
+        lines_by_intensity = df1.S.tolist()
+        lines_by_intensity.sort()
 
         if len(df) == 0:  # no lines
             self._Nlines_cutoff = None
@@ -3212,6 +3222,53 @@ def _cutoff_linestrength(self, cutoff=None):
 
             error = df.S[b].sum() / df.S.sum() * 100
 
+            if cutoff_error is not None:
+                if cutoff_error < error:
+                    # Remove additional lines such that cutoff error is less than user-inputted value
+                    intensity_sum = 0
+                    for i in lines_by_intensity:
+                        intensity_sum += i
+                    # Last value of the list is the cumulative intensity
+                    lines_by_intensity.append(intensity_sum)
+                    desired_intensity = (cutoff_error / 100) * lines_by_intensity[-1]
+                    next_index = 0
+                    sum = 0
+                    # Find the index till which sum <= desired intensity and remove the rest of the elements
+                    for i in lines_by_intensity[:-1]:
+                        if sum <= desired_intensity:
+                            if sum + lines_by_intensity[next_index] > desired_intensity:
+                                break
+                            else:
+                                sum += i
+                                next_index += 1
+                        else:
+                            break
+                    # Change a copy of the list containing lines by intensity
+                    correct_lines = lines_by_intensity.copy()
+                    for i in range(next_index + 1, len(lines_by_intensity) + 1):
+                        correct_lines.pop()
+                    # Update the database
+                    df1["S"] = correct_lines
+                    # Print current error
+                    current_error = sum / lines_by_intensity[-1] * 100
+                    print(
+                        "Current percentage error: {0:.2f}% ".format(current_error)
+                        + "Inputted error: {0:.2f}%".format(cutoff_error)
+                    )
+
+                else:
+                    print(
+                        "Cutoff error is greater than the current error."
+                        + " Inputted cutoff error percentage is {0:.2f}%".format(
+                            cutoff_error
+                        )
+                        + " Estimated error: {0:.2f}%".format(error)
+                    )
+                b = df.S <= cutoff
+                Nlines_cutoff = b.sum()
+            else:
+                print("Cutoff error not inputted.")
+
             if verbose >= 2:
                 print(
                     "Discarded {0:.2f}% of lines (linestrength<{1}cm-1/(#.cm-2))".format(

radis/lbl/loader.py

@@ -451,6 +451,7 @@ class Parameters(ConditionDict):
         "neighbour_lines",
         "chunksize",
         "cutoff",
+        "cutoff_error",
         "db_use_cached",
         "dbformat",
         "dbpath",
@@ -483,6 +484,9 @@ def __init__(self):
         self.truncation = None  #: float: cutoff for half-width lineshape calculation (cm-1). Overwritten by SpectrumFactory
         self.neighbour_lines = None  #: float: extra range (cm-1) on each side of the spectrum to account for neighbouring lines. Overwritten by SpectrumFactory
         self.cutoff = None  #: float: linestrength cutoff (molecule/cm)
+        self.cutoff_error = (
+            None  #: float: error percentage below which the cutoff error should be
+        )
         self.broadening_method = ""  #: str:``"voigt"``, ``"convolve"``, ``"fft"``
         self.optimization = None  #: str: ``"simple"``, ``"min-RMS"``, ``None``
         self.db_use_cached = (