Taking the better range finding functions from BP and DB analyses, added a size check.

Only affects one or two pixels per module.

CS_analysis.cpp

@@ -18,36 +18,52 @@
 #include <algorithm>
 #include <sstream>
 
-double lowestFullySwitchedPoint(const vector<double> &higher_filter, const vector<double> &lower_filter) {
+double highestPointBeforeSwitching(const vector<double> &lower_filter, const vector<double> &higher_filter) {
 
-  // find the lowest value in higher_filter that does not appear in lower_filter
+  double highest_point_before_switching;
-  double lowest_fully_switched_val = *max_element(higher_filter.begin(),higher_filter.end())+1;
 
-  for(vector<double>::const_iterator it = higher_filter.begin(); it != higher_filter.end(); ++it) {
+  if (higher_filter.size() > 0) {
-    if (find(lower_filter.begin(), lower_filter.end(), *it) == lower_filter.end()) {
+    // find the highest value in lower_filter that is lower than all entries in higher_filter
-      if (*it < lowest_fully_switched_val) {
+    highest_point_before_switching = *min_element(lower_filter.begin(),lower_filter.end())-1;
-	lowest_fully_switched_val = *it;
-      }
-    }
-  }
-
-  return lowest_fully_switched_val;
-}
-
-double highestFullySwitchedPoint(const vector<double> &lower_filter, const vector<double> &higher_filter) {
-
-  // find the highest value in lower_filter that does not appear in higher_filter
-  double highest_fully_switched_val = *min_element(lower_filter.begin(),lower_filter.end())-1;
     
+    double lowest_entry_in_higher_filter = *min_element(higher_filter.begin(),higher_filter.end());
     for(vector<double>::const_iterator it = lower_filter.begin(); it != lower_filter.end(); ++it) {
-    if (find(higher_filter.begin(), higher_filter.end(), *it) == higher_filter.end()) {
+      if (*it < lowest_entry_in_higher_filter) {
-      if (*it > highest_fully_switched_val) {
+	if (*it > highest_point_before_switching) {
-	highest_fully_switched_val = *it;
+	  highest_point_before_switching = *it;
 	}
       }
     }
+  } else {
+    highest_point_before_switching = *max_element(lower_filter.begin(),lower_filter.end());
+  }
+
+  return highest_point_before_switching;
+
+}
+
+double lowestPointAfterSwitching(const vector<double> &higher_filter, const vector<double> &lower_filter) {
+
+  double lowest_point_after_switching;
+
+  if (lower_filter.size() > 0) {
+    // find the lowest value in higher_filter that is higher than all entries in lower_filter
+    lowest_point_after_switching = *max_element(higher_filter.begin(),higher_filter.end())+1;
+
+    double highest_entry_in_lower_filter = *max_element(lower_filter.begin(),lower_filter.end());
+    for(vector<double>::const_iterator it = higher_filter.begin(); it != higher_filter.end(); ++it) {
+      if (*it > highest_entry_in_lower_filter) {
+	if (*it < lowest_point_after_switching) {
+	  lowest_point_after_switching = *it;
+	}
+      }
+    }
+  } else {
+    lowest_point_after_switching = *min_element(higher_filter.begin(),higher_filter.end());
+  }
+
+  return lowest_point_after_switching;
 
-  return highest_fully_switched_val;
 }
 
 double checkRangeMinForRCEffect(double range_min) {
@@ -1081,7 +1097,7 @@ int main(int argc, char* argv[]) {
 	rangemin_g0 = *min_element(r0_filter.begin(),r0_filter.end());
 	rangemin_g0 = checkRangeMinForRCEffect(rangemin_g0);
 
-	rangemax_g0 = highestFullySwitchedPoint(r0_filter, r1_filter);
+	rangemax_g0 = highestPointBeforeSwitching(r0_filter, r1_filter);
 
 	if (rangemax_g0 > rangemin_g0) {
 
@@ -1173,10 +1189,10 @@ int main(int argc, char* argv[]) {
 
       if (r1_adc.size() > 1) {
 
-	rangemin_g1 = lowestFullySwitchedPoint(r1_filter, r0_filter);
+	rangemin_g1 = lowestPointAfterSwitching(r1_filter, r0_filter);
 	rangemin_g1 = checkRangeMinForRCEffect(rangemin_g1);
 
-	rangemax_g1 = highestFullySwitchedPoint(r1_filter, r2_filter);
+	rangemax_g1 = highestPointBeforeSwitching(r1_filter, r2_filter);
 
 	if (rangemax_g1 > rangemin_g1) {
 
@@ -1293,7 +1309,7 @@ int main(int argc, char* argv[]) {
 
       if (r2_adc.size() > 1) {
 
-	rangemin_g2 = lowestFullySwitchedPoint(r2_filter, r1_filter);
+	rangemin_g2 = lowestPointAfterSwitching(r2_filter, r1_filter);
 	rangemin_g2 = checkRangeMinForRCEffect(rangemin_g2);
 
       	rangemax_g2 = *max_element(r2_filter.begin(),r2_filter.end());
@@ -1553,6 +1569,7 @@ int main(int argc, char* argv[]) {
     g0overg1hist->Draw();
     c1->Update();
     TF1* g0overg1hist_gaus = g0overg1hist->GetFunction("gaus");
+    if (g0overg1hist_gaus) {
       sprintf(savename,"Mean %0.4f #pm %0.4f", g0overg1hist_gaus->GetParameter(1), g0overg1hist_gaus->GetParError(1));
       pave2->AddText((const char *)(savename));
       sprintf(savename,"Sigma %0.4f #pm %0.4f", g0overg1hist_gaus->GetParameter(2), g0overg1hist_gaus->GetParError(2));
@@ -1561,6 +1578,7 @@ int main(int argc, char* argv[]) {
       g0overg1hist->SetStats(kFALSE);
       sprintf(savename,"plots/M%s/CurrentSource/g0overg1hist_fit_M%s.png", module_str.c_str(), module_str.c_str());
       c1->SaveAs((const char *)(savename));
+    }
 
     g0overg1erhist->GetXaxis()->SetTitle("G0 / G1 uncert");
     g0overg1erhist->GetXaxis()->SetTitleOffset(1.1);
@@ -1621,6 +1639,7 @@ int main(int argc, char* argv[]) {
     g1overg2hist->Draw();
     c1->Update();
     TF1* g1overg2hist_gaus = g1overg2hist->GetFunction("gaus");
+    if (g1overg2hist_gaus) {
       sprintf(savename,"Mean %0.4f #pm %0.4f", g1overg2hist_gaus->GetParameter(1), g1overg2hist_gaus->GetParError(1));
       pave2->AddText((const char *)(savename));
       sprintf(savename,"Sigma %0.4f #pm %0.4f", g1overg2hist_gaus->GetParameter(2), g1overg2hist_gaus->GetParError(2));
@@ -1629,6 +1648,7 @@ int main(int argc, char* argv[]) {
       g1overg2hist->SetStats(kFALSE);
       sprintf(savename,"plots/M%s/CurrentSource/g1overg2hist_fit_M%s.png", module_str.c_str(), module_str.c_str());
       c1->SaveAs((const char *)(savename));
+    }
 
     g1overg2erhist->GetXaxis()->SetTitle("G1 / G2 uncert");
     g1overg2erhist->GetXaxis()->SetTitleOffset(1.1);