update GaussHitFinderSBN_module to take per plane hit finding parameters

Author: Michael Carrigan

sbncode/CAFMaker/RecoUtils/CMakeLists.txt

@@ -11,4 +11,6 @@ art_make_library( LIBRARY_NAME caf_RecoUtils
                     larsim::MCCheater_BackTrackerService_service
                     larsim::MCCheater_ParticleInventoryService_service
                     larcorealg::Geometry
+                    sbnanaobj::StandardRecord
+                    sbnanaobj::StandardRecordFlat
                 )

sbncode/HitFinder/GaussHitFinderSBN_module.cc

@@ -79,11 +79,11 @@ namespace hit {
     const std::vector<int> fLongMaxHitsVec;    ///<Maximum number hits on a really long pulse train
     const std::vector<int> fLongPulseWidthVec; ///<Sets width of hits used to describe long pulses
 
-    const size_t fMaxMultiHit; ///<maximum hits for multi fit
+    const std::vector<size_t> fMaxMultiHit; ///<maximum hits for multi fit
     const int fAreaMethod;     ///<Type of area calculation
     const std::vector<double>
       fAreaNormsVec;       ///<factors for converting area to same units as peak height
-    const double fChi2NDF; ///maximum Chisquared / NDF allowed for a hit to be saved
+    const std::vector<double> fChi2NDF; ///maximum Chisquared / NDF allowed for a hit to be saved
 
     const std::vector<float> fPulseHeightCuts;
     const std::vector<float> fPulseWidthCuts;
@@ -116,10 +116,10 @@ namespace hit {
     , fLongMaxHitsVec(pset.get<std::vector<int>>("LongMaxHits", std::vector<int>() = {25, 25, 25}))
     , fLongPulseWidthVec(
         pset.get<std::vector<int>>("LongPulseWidth", std::vector<int>() = {16, 16, 16}))
-    , fMaxMultiHit(pset.get<int>("MaxMultiHit"))
+    , fMaxMultiHit(pset.get<std::vector<size_t>>("MaxMultiHitPerPlane", std::vector<size_t>() = {5, 5, 5}))
     , fAreaMethod(pset.get<int>("AreaMethod"))
     , fAreaNormsVec(FillOutHitParameterVector(pset.get<std::vector<double>>("AreaNorms")))
-    , fChi2NDF(pset.get<double>("Chi2NDF"))
+    , fChi2NDF(pset.get<std::vector<double>>("Chi2NDFPerPlane", std::vector<double>() = {500.0, 500.0, 500.0}))
     , fPulseHeightCuts(
         pset.get<std::vector<float>>("PulseHeightCuts", std::vector<float>() = {3.0, 3.0, 3.0}))
     , fPulseWidthCuts(
@@ -365,13 +365,13 @@ namespace hit {
               // #######################################################
               // ### If # requested Gaussians is too large then punt ###
               // #######################################################
-              if (mergedCands.size() <= fMaxMultiHit) {
+              if (mergedCands.size() <= fMaxMultiHit.at(plane)) {
                 fPeakFitterTool->findPeakParameters(
                   range.data(), mergedCands, peakParamsVec, chi2PerNDF, NDF);
 
                 // If the chi2 is infinite then there is a real problem so we bail
                 if (!(chi2PerNDF < std::numeric_limits<double>::infinity())) {
-                  chi2PerNDF = 2. * fChi2NDF;
+                  chi2PerNDF = 2. * fChi2NDF.at(plane);
                   NDF = 2;
                 }
 
@@ -384,7 +384,7 @@ namespace hit {
               // ###   depend on the fhicl parameter fLongPulseWidth ###
               // ### Also do this if chi^2 is too large              ###
               // #######################################################
-              if (mergedCands.size() > fMaxMultiHit || nGausForFit * chi2PerNDF > fChi2NDF) {
+              if (mergedCands.size() > fMaxMultiHit.at(plane) || nGausForFit * chi2PerNDF > fChi2NDF.at(plane)) {
                 int longPulseWidth = fLongPulseWidthVec.at(plane);
                 int nHitsThisPulse = (endT - startT) / longPulseWidth;
 
@@ -401,7 +401,7 @@ namespace hit {
                 peakParamsVec.clear();
                 nGausForFit = nHitsThisPulse;
                 NDF = 1.;
-                chi2PerNDF = chi2PerNDF > fChi2NDF ? chi2PerNDF : -1.;
+                chi2PerNDF = chi2PerNDF > fChi2NDF.at(plane) ? chi2PerNDF : -1.;
 
                 for (int hitIdx = 0; hitIdx < nHitsThisPulse; hitIdx++) {
                   // This hit parameters

sbncode/HitFinder/hitfindermodules_sbn.fcl

@@ -11,17 +11,19 @@ gauss_hitfinder:
 {
     module_type:          "GaussHitFinderSBN"
     CalDataModuleLabel:   "caldata"
-    MaxMultiHit:          5                  # maximum hits for multi gaussia fit attempt
-    AreaMethod:           0                  # 0 = area by integral, 1 = area by gaussian area formula
-    AreaNorms:            [ 1.0, 1.0, 1.0 ]  # normalizations that put signal area in
-                                             # same scale as peak height.
-    TryNplus1Fits:        false              # Don't try to refit with extra peak if bad chisq
-    LongMaxHits:          [ 25, 25, 25]      # max number hits in long pulse trains
-    LongPulseWidth:       [ 10, 10, 10]      # max widths for hits in long pulse trains
-    Chi2NDF:              500.               # maximum Chisquared / NDF allowed to store fit, if fail
-                                             # will use "long" pulse method to return hit
-    AllHitsInstanceName:  ""                 # If non-null then this will be the instance name of all hits output to event
-                                             # in this case there will be two hit collections, one filtered and one containing all hits
+    MaxMultiHit:          5                       # maximum hits for multi gaussian fit attempt
+    MaxMultiHitPerPlane:  [ 10, 10, 8 ]           # maximum hits per plane for multi gaussia fit attempt
+    AreaMethod:           0                       # 0 = area by integral, 1 = area by gaussian area formula
+    AreaNorms:            [ 1.0, 1.0, 1.0 ]       # normalizations that put signal area in
+                                                  # same scale as peak height.
+    TryNplus1Fits:        false                   # Don't try to refit with extra peak if bad chisq
+    LongMaxHits:          [ 25, 25, 25]           # max number hits in long pulse trains
+    LongPulseWidth:       [ 5, 5, 5]              # max widths for hits in long pulse trains
+    Chi2NDF:              500.                    # maximum Chisquared / NDF allowed to store fit
+    Chi2NDFPerPlane:      [ 2500., 1750., 2500.]  # maximum Chisquared / NDF allowed per plane to store, if fail
+                                                  # will use "long" pulse method to return hit
+    AllHitsInstanceName:  ""                      # If non-null then this will be the instance name of all hits output to event
+                                                  # in this case there will be two hit collections, one filtered and one containing all hits
 
     # Candididate peak finding done by tool, one tool instantiated per plane (but could be other divisions too)
     HitFinderToolVec:
@@ -39,7 +41,7 @@ gauss_hitfinder:
     HitFilterAlg:
     {
         AlgName: "HitFilterAlg"
-        MinPulseHeight: [5.0, 5.0, 5.0]      #minimum hit peak amplitude per plane
+        MinPulseHeight: [2.0, 2.0, 2.0]      #minimum hit peak amplitude per plane
         MinPulseSigma:  [1.0, 1.0, 1.0]      #minimum hit rms per plane
     }
                                              # In addition to the filter alg we can also filter hits on the same pulse train