Skip to content

Commit 2e800b2

Browse files
committed
Merge branch 'feature/carriganm95_hitTuning' into develop
2 parents bd52c28 + cf2f906 commit 2e800b2

2 files changed

Lines changed: 43 additions & 24 deletions

File tree

sbncode/HitFinder/GaussHitFinderSBN_module.cc

Lines changed: 31 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -69,21 +69,39 @@ namespace hit {
6969
void beginJob(art::ProcessingFrame const&) override;
7070

7171
std::vector<double> FillOutHitParameterVector(const std::vector<double>& input);
72+
73+
template<class T>
74+
std::vector<T> getValueOrListOf(fhicl::ParameterSet const& pset, std::string const& key) const {
75+
76+
auto const& wireReadoutGeom = art::ServiceHandle<geo::WireReadout const>()->Get();
77+
const unsigned int N_PLANES = wireReadoutGeom.Nplanes();
78+
79+
if (pset.is_key_to_sequence(key)) {
80+
auto vec = pset.get<std::vector<T>>(key);
81+
if (vec.size() != N_PLANES) {
82+
throw cet::exception("Configuration")
83+
<< "Parameter '" << key << "' has size " << vec.size()
84+
<< " but expected " << N_PLANES << " (number of planes)";
85+
}
86+
return pset.get<std::vector<T>>(key);
87+
} else
88+
return std::vector<T>(N_PLANES, pset.get<T>(key));
89+
} // getValueOrListOf()
7290

7391
const bool fFilterHits;
7492
const bool fFillHists;
7593

7694
const std::string fCalDataModuleLabel;
7795
const std::string fAllHitsInstanceName;
7896

79-
const std::vector<int> fLongMaxHitsVec; ///<Maximum number hits on a really long pulse train
80-
const std::vector<int> fLongPulseWidthVec; ///<Sets width of hits used to describe long pulses
97+
const std::vector<int> fLongMaxHitsVec; ///< Maximum number hits on a really long pulse train
98+
const std::vector<int> fLongPulseWidthVec; ///< Sets width of hits used to describe long pulses
8199

82-
const size_t fMaxMultiHit; ///<maximum hits for multi fit
83-
const int fAreaMethod; ///<Type of area calculation
100+
const std::vector<size_t> fMaxMultiHit; ///< Maximum hits for multi fit
101+
const int fAreaMethod; ///< Type of area calculation
84102
const std::vector<double>
85-
fAreaNormsVec; ///<factors for converting area to same units as peak height
86-
const double fChi2NDF; ///maximum Chisquared / NDF allowed for a hit to be saved
103+
fAreaNormsVec; ///< Factors for converting area to same units as peak height
104+
const std::vector<double> fChi2NDF; ///< Maximum Chisquared / NDF allowed for a hit to be saved
87105

88106
const std::vector<float> fPulseHeightCuts;
89107
const std::vector<float> fPulseWidthCuts;
@@ -116,17 +134,18 @@ namespace hit {
116134
, fLongMaxHitsVec(pset.get<std::vector<int>>("LongMaxHits", std::vector<int>() = {25, 25, 25}))
117135
, fLongPulseWidthVec(
118136
pset.get<std::vector<int>>("LongPulseWidth", std::vector<int>() = {16, 16, 16}))
119-
, fMaxMultiHit(pset.get<int>("MaxMultiHit"))
137+
, fMaxMultiHit(getValueOrListOf<size_t>(pset, "MaxMultiHit"))
120138
, fAreaMethod(pset.get<int>("AreaMethod"))
121139
, fAreaNormsVec(FillOutHitParameterVector(pset.get<std::vector<double>>("AreaNorms")))
122-
, fChi2NDF(pset.get<double>("Chi2NDF"))
140+
, fChi2NDF(getValueOrListOf<double>(pset, "Chi2NDF"))
123141
, fPulseHeightCuts(
124142
pset.get<std::vector<float>>("PulseHeightCuts", std::vector<float>() = {3.0, 3.0, 3.0}))
125143
, fPulseWidthCuts(
126144
pset.get<std::vector<float>>("PulseWidthCuts", std::vector<float>() = {2.0, 1.5, 1.0}))
127145
, fPulseRatioCuts(
128146
pset.get<std::vector<float>>("PulseRatioCuts", std::vector<float>() = {0.35, 0.40, 0.20}))
129147
{
148+
130149
if (fFillHists && art::Globals::instance()->nthreads() > 1u) {
131150
throw art::Exception(art::errors::Configuration)
132151
<< "Cannot fill histograms when multiple threads configured, please set fFillHists to "
@@ -365,13 +384,13 @@ namespace hit {
365384
// #######################################################
366385
// ### If # requested Gaussians is too large then punt ###
367386
// #######################################################
368-
if (mergedCands.size() <= fMaxMultiHit) {
387+
if (mergedCands.size() <= fMaxMultiHit.at(plane)) {
369388
fPeakFitterTool->findPeakParameters(
370389
range.data(), mergedCands, peakParamsVec, chi2PerNDF, NDF);
371390

372391
// If the chi2 is infinite then there is a real problem so we bail
373392
if (!(chi2PerNDF < std::numeric_limits<double>::infinity())) {
374-
chi2PerNDF = 2. * fChi2NDF;
393+
chi2PerNDF = 2. * fChi2NDF.at(plane);
375394
NDF = 2;
376395
}
377396

@@ -384,7 +403,7 @@ namespace hit {
384403
// ### depend on the fhicl parameter fLongPulseWidth ###
385404
// ### Also do this if chi^2 is too large ###
386405
// #######################################################
387-
if (mergedCands.size() > fMaxMultiHit || nGausForFit * chi2PerNDF > fChi2NDF) {
406+
if (mergedCands.size() > fMaxMultiHit.at(plane) || nGausForFit * chi2PerNDF > fChi2NDF.at(plane)) {
388407
int longPulseWidth = fLongPulseWidthVec.at(plane);
389408
int nHitsThisPulse = (endT - startT) / longPulseWidth;
390409

@@ -401,7 +420,7 @@ namespace hit {
401420
peakParamsVec.clear();
402421
nGausForFit = nHitsThisPulse;
403422
NDF = 1.;
404-
chi2PerNDF = chi2PerNDF > fChi2NDF ? chi2PerNDF : -1.;
423+
chi2PerNDF = chi2PerNDF > fChi2NDF.at(plane) ? chi2PerNDF : -1.;
405424

406425
for (int hitIdx = 0; hitIdx < nHitsThisPulse; hitIdx++) {
407426
// This hit parameters

sbncode/HitFinder/hitfindermodules_sbn.fcl

Lines changed: 12 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -11,17 +11,17 @@ gauss_hitfinder:
1111
{
1212
module_type: "GaussHitFinderSBN"
1313
CalDataModuleLabel: "caldata"
14-
MaxMultiHit: 5 # maximum hits for multi gaussia fit attempt
15-
AreaMethod: 0 # 0 = area by integral, 1 = area by gaussian area formula
16-
AreaNorms: [ 1.0, 1.0, 1.0 ] # normalizations that put signal area in
17-
# same scale as peak height.
18-
TryNplus1Fits: false # Don't try to refit with extra peak if bad chisq
19-
LongMaxHits: [ 25, 25, 25] # max number hits in long pulse trains
20-
LongPulseWidth: [ 10, 10, 10] # max widths for hits in long pulse trains
21-
Chi2NDF: 500. # maximum Chisquared / NDF allowed to store fit, if fail
22-
# will use "long" pulse method to return hit
23-
AllHitsInstanceName: "" # If non-null then this will be the instance name of all hits output to event
24-
# in this case there will be two hit collections, one filtered and one containing all hits
14+
MaxMultiHit: 5 # maximum hits for multi gaussian fit attempt
15+
AreaMethod: 0 # 0 = area by integral, 1 = area by gaussian area formula
16+
AreaNorms: [ 1.0, 1.0, 1.0 ] # normalizations that put signal area in
17+
# same scale as peak height.
18+
TryNplus1Fits: false # Don't try to refit with extra peak if bad chisq
19+
LongMaxHits: [ 25, 25, 25] # max number hits in long pulse trains
20+
LongPulseWidth: [ 5, 5, 5] # max widths for hits in long pulse trains
21+
Chi2NDF: 500. # maximum Chisquared / NDF allowed to store fit
22+
# will use "long" pulse method to return hit
23+
AllHitsInstanceName: "" # If non-null then this will be the instance name of all hits output to event
24+
# in this case there will be two hit collections, one filtered and one containing all hits
2525

2626
# Candididate peak finding done by tool, one tool instantiated per plane (but could be other divisions too)
2727
HitFinderToolVec:
@@ -39,7 +39,7 @@ gauss_hitfinder:
3939
HitFilterAlg:
4040
{
4141
AlgName: "HitFilterAlg"
42-
MinPulseHeight: [5.0, 5.0, 5.0] #minimum hit peak amplitude per plane
42+
MinPulseHeight: [2.0, 2.0, 2.0] #minimum hit peak amplitude per plane
4343
MinPulseSigma: [1.0, 1.0, 1.0] #minimum hit rms per plane
4444
}
4545
# In addition to the filter alg we can also filter hits on the same pulse train

0 commit comments

Comments
 (0)