fix: "WARNING in HistoUtil.getHistoIQR: attempt to calculate median o…

…f an empty list" warnings from RG-A Sp18 (#225)
JeffersonLab · Aug 19, 2024 · fcbb500 · fcbb500
1 parent 4f0ac7d
commit fcbb500
Show file tree

Hide file tree

Showing 4 changed files with 19 additions and 9 deletions.
diff --git a/bin/run-detectors-timelines.sh b/bin/run-detectors-timelines.sh
@@ -338,7 +338,7 @@ else
 fi
 
 # grep for suspicious things in error logs
-errPattern="error:|exception:"
+errPattern="error:|exception:|warning"
 echo """To look for any quieter errors, running \`grep -iE '$errPattern'\` on *.err files:
 $sep"""
 grep -iE --color "$errPattern" $logDir/*.err || echo "Good news: grep found no errors, but you still may want to take a look yourself..."

diff --git a/detectors/src/main/java/org/jlab/clas/timeline/fitter/CTOFFitter.groovy b/detectors/src/main/java/org/jlab/clas/timeline/fitter/CTOFFitter.groovy
@@ -10,7 +10,7 @@ class CTOFFitter {
     def f1 = new F1D("fit:"+h1.getName(), "[amp]*gaus(x,[mean],[sigma])", -5.0, 5.0)
     double hAmp   = h1.getBinContent(h1.getMaximumBin());
     double hMu    = h1.getAxis().getBinCenter(h1.getMaximumBin())
-    double hSigma = Math.min(HistoUtil.getHistoIQR(h1) / 2.0, 0.2)
+    double hSigma = Math.min(HistoUtil.getHistoIQR(h1, 0.4) / 2.0, 0.2)
     def rangeFactor = 1.5
     double rangeMin = (hMu - rangeFactor*hSigma);
     double rangeMax = (hMu + rangeFactor*hSigma);

diff --git a/detectors/src/main/java/org/jlab/clas/timeline/fitter/HTCCFitter.groovy b/detectors/src/main/java/org/jlab/clas/timeline/fitter/HTCCFitter.groovy
@@ -18,7 +18,7 @@ class HTCCFitter{
 
     // compute the IQR of `h1`, which is a better estimate of the initial sigma parameter than
     // the RMS, since it is not (typically) biased by outliers
-    def iqr = HistoUtil.getHistoIQR(h1)
+    def iqr = HistoUtil.getHistoIQR(h1, -1.0)
 
     def f1 = new F1D("fit:"+h1.getName(), "[amp]*gaus(x,[mean],[sigma])", maxV-2, maxV+3);
 
@@ -28,8 +28,13 @@ class HTCCFitter{
 
     f1.setParameter(0, h1.getMax());
     f1.setParameter(1, maxV);
-    f1.setParameter(2, Math.max(iqr / 2.0, 0.1));
-    f1.setParLimits(2, Math.max(iqr / 5.0, 0.05), Math.min(iqr * 5.0, h1.getAxis().getBinCenter(h1.getAxis().getNBins()-1)));
+    if(iqr > 0) {
+      f1.setParameter(2, Math.max(iqr / 2.0, 0.1));
+      f1.setParLimits(2, Math.max(iqr / 5.0, 0.05), Math.min(iqr * 5.0, h1.getAxis().getBinCenter(h1.getAxis().getNBins()-1)));
+    } else {
+      f1.setParameter(2, 0.1);
+      f1.setParLimits(2, 0.05, h1.getAxis().getBinCenter(h1.getAxis().getNBins()-1));
+    }
 
     DataFitter.fit(f1, h1, "");
 

diff --git a/detectors/src/main/java/org/jlab/clas/timeline/util/HistoUtil.groovy b/detectors/src/main/java/org/jlab/clas/timeline/util/HistoUtil.groovy
@@ -59,8 +59,12 @@ class HistoUtil {
   }
 
   /// @returns the Interquartile Range (IQR) of a histogram
-  /// @param the histogram; FIXME: only unweighted histograms are supported at the moment
-  static def getHistoIQR(H1F hist) {
+  /// @param hist the histogram; FIXME: only unweighted histograms are supported at the moment
+  /// @param defIfEmpty the default value, if `hist` is empty
+  static def getHistoIQR(H1F hist, defIfEmpty) {
+    if( !(hist.getEntries() > 0) || !(hist.integral() > 0) ) { // `.getEntries()` can be nonzero for empty histograms, so be sure to check the integral too
+      return defIfEmpty
+    }
     def hist_list = []
     hist.getAxis().getNBins().times { bin ->
       def counts = hist.getBinContent(bin).toInteger() // FIXME: assumes the histogram is unweighted
@@ -69,8 +73,9 @@ class HistoUtil {
     }
     def listMedian = { d ->
       if(d.size()==0) {
-        System.err.println("WARNING in HistoUtil.getHistoIQR: attempt to calculate median of an empty list")
-        return 0
+        // this list may end up being empty if there are _few_ entries in the histogram, since this method
+        // is called to get the quartiles; in this case, just return `defIfEmpty`
+        return defIfEmpty
       }
       d.sort()
       def m = d.size().intdiv(2)