Belle II Software development
TrackIsoCalculatorModule.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#include <analysis/modules/TrackIsoCalculator/TrackIsoCalculatorModule.h>
10
11#include <analysis/dbobjects/PIDDetectorWeights.h>
12#include <analysis/DecayDescriptor/DecayDescriptorParticle.h>
13#include <analysis/VariableManager/Manager.h>
14#include <analysis/utility/DetectorSurface.h>
15
16#include <framework/database/DBObjPtr.h>
17
18#include <cmath>
19#include <boost/algorithm/string.hpp>
20
21
22using namespace Belle2;
23
24REG_MODULE(TrackIsoCalculator);
25
27{
28 // Set module properties
30 R"DOC(Calculate track isolation variables on the charged stable particles, or selected charged daughters, of the input ParticleList.)DOC");
32
33 // Parameter definitions
34 addParam("decayString",
36 "The name of the input charged stable particle list, e.g. ``mu+:all``, or a composite particle w/ charged stable daughters for which distances are to be calculated, e.g. ``D0 -> ^K- pi+``. Note that in the latter case we allow only one daughter to be selected in the decay string per module instance.");
37 addParam("particleListReference",
39 "The name of the input ParticleList of reference tracks. Must be a charged stable particle as defined in Const::chargedStableSet.");
40 addParam("detectorNames",
42 "The list of names of the detectors at whose (cylindrical) surface(s) we extrapolate each helix's polar and azimuthal angle. Allowed values: {CDC, TOP, ARICH, ECL, KLM}.",
43 std::vector<std::string>());
44 addParam("useHighestProbMassForExt",
46 "If this option is set, the helix extrapolation for the target and reference particles will use the track fit result for the most probable mass hypothesis, namely, the one that gives the highest chi2Prob of the fit.",
47 bool(false));
48 addParam("payloadName",
50 "The name of the database payload object with the PID detector weights.",
51 std::string("PIDDetectorWeights"));
52 addParam("excludePIDDetWeights",
54 "If set to true, will not use the PID detector weights for the score definition.",
55 bool(false));
56}
57
61
63{
64 m_event_metadata.isRequired();
65
67 m_DBWeights = std::make_unique<DBObjPtr<PIDDetectorWeights>>(m_payloadName);
68 }
69
70 bool valid = m_decaydescriptor.init(m_decayString);
71 if (!valid) {
72 B2ERROR("Decay string " << m_decayString << " is invalid.");
73 }
74
75 const DecayDescriptorParticle* mother = m_decaydescriptor.getMother();
76
77 m_pListTarget.isRequired(mother->getFullName());
78 m_nSelectedDaughters = m_decaydescriptor.getSelectionNames().size();
79
80 if (m_nSelectedDaughters > 1) {
81 B2ERROR("More than one daughter is selected in the decay string " << m_decayString << ".");
82 }
83
85
86 B2INFO("Calculating track-based isolation variables for the decay string: "
88 << " using the reference ParticleList: "
90 << ".");
91
93 B2ERROR("Selected ParticleList in decay string:"
95 << " and/or reference ParticleList: "
97 << " is not that of a valid particle in Const::chargedStableSet!");
98 }
99
101 B2INFO("Will use track fit result for the most probable mass hypothesis in helix extrapolation.");
102 }
103
104 std::string detNamesConcat("");
105
106 for (auto& detName : m_detNames) {
107
108 // Ensure we use uppercase detector labels.
109 boost::to_upper(detName);
110
111 detNamesConcat += "_" + detName;
112
113 // Define the name(s) of the variables for this detector to be stored as extraInfo.
114 for (const auto& iLayer : DetectorSurface::detToLayers.at(detName)) {
115
116
117 auto iDetLayer = detName + std::to_string(iLayer);
118
119 auto distVarName = "distToClosestTrkAt" + iDetLayer + "_VS_" + m_pListReferenceName;
121 distVarName += "__useHighestProbMassForExt";
122 }
123 auto refPartIdxVarName = "idxOfClosestPartAt" + iDetLayer + "In_" + m_pListReferenceName;
124
125 m_detLayerToDistVariable.insert(std::make_pair(iDetLayer, distVarName));
126 m_detLayerToRefPartIdxVariable.insert(std::make_pair(iDetLayer, refPartIdxVarName));
127
128 }
129
130 // Isolation score variable.
131 m_isoScoreVariable = "trkIsoScore" + detNamesConcat + "_VS_" + m_pListReferenceName;
132 m_isoScoreVariableAsWeightedAvg = "trkIsoScoreAsWeightedAvg" + detNamesConcat + "_VS_" + m_pListReferenceName;
134 m_isoScoreVariable += "__useHighestProbMassForExt";
135 m_isoScoreVariableAsWeightedAvg += "__useHighestProbMassForExt";
136 }
137 }
138
139}
140
142{
143 B2DEBUG(11, "Start processing EVENT: " << m_event_metadata->getEvent());
144
145 const auto nMotherParticles = m_pListTarget->getListSize();
146
147 // Fill transient container of all selected charged particles in the decay
148 // for which the distance to nearest neighbour is to be calculated.
149 // If the input ParticleList is that of standard charged particles, just copy
150 // the whole list over, otherwise loop over the selected charged daughters.
151 std::unordered_map<int, const Particle*> targetParticles;
152 targetParticles.reserve(nMotherParticles);
153
154 for (unsigned int iPart(0); iPart < nMotherParticles; ++iPart) {
155
156 auto iParticle = m_pListTarget->getParticle(iPart);
157
158 // Loop over the input detectors.
159 for (const auto& detName : m_detNames) {
160
161 // Loop over the layer of the input detector.
162 for (const auto& iLayer : DetectorSurface::detToLayers.at(detName)) {
163
164 auto iDetLayer = detName + std::to_string(iLayer);
165
167 for (auto* iDaughter : m_decaydescriptor.getSelectionParticles(iParticle)) {
168 // Check if the distance for this target particle has been set already,
169 // e.g. by a previous instance of this module.
170 if (iDaughter->hasExtraInfo(m_detLayerToDistVariable[iDetLayer])) {
171 continue;
172 }
173 targetParticles.insert({iDaughter->getMdstSource(), iDaughter});
174 }
175 } else {
176 if (iParticle->hasExtraInfo(m_detLayerToDistVariable[iDetLayer])) {
177 continue;
178 }
179 targetParticles.insert({iParticle->getMdstSource(), iParticle});
180 }
181
182 }
183 }
184 }
185
186 const auto nParticlesTarget = targetParticles.size();
187 const auto nParticlesReference = m_pListReference->getListSize();
188
189 // Loop over the input detectors.
190 for (const auto& detName : m_detNames) {
191
192 // Loop over the layer of the input detector.
193 for (const auto& iLayer : DetectorSurface::detToLayers.at(detName)) {
194
195 auto iDetLayer = detName + std::to_string(iLayer);
196
197 B2DEBUG(11, "\n"
198 << "Detector surface: " << iDetLayer << "\n"
199 << "nMotherParticles: " << nMotherParticles << "\n"
200 << "nParticlesTarget: " << nParticlesTarget << "\n"
201 << "nParticlesReference: " << nParticlesReference);
202
203 double dummyDist(-1.0);
204
205 // Store the pair-wise distances in a map,
206 // where the keys are pairs of mdst indexes.
207 std::map<std::pair<int, int>, double> particleMdstSourcePairsToDist;
208
209 // Loop over input particle list
210 for (const auto& targetParticle : targetParticles) {
211
212 auto iMdstSource = targetParticle.first;
213 auto iParticle = targetParticle.second;
214
215 for (unsigned int jPart(0); jPart < nParticlesReference; ++jPart) {
216
217 auto jParticle = m_pListReference->getParticle(jPart);
218 auto jMdstSource = jParticle->getMdstSource();
219
220 auto partMdstSourcePair = std::make_pair(iMdstSource, jMdstSource);
221
222 // Set dummy distance if same particle.
223 if (iMdstSource == jMdstSource) {
224 particleMdstSourcePairsToDist[partMdstSourcePair] = dummyDist;
225 continue;
226 }
227 // If:
228 //
229 // - the mass hypothesis of the best fit is used, OR
230 // - the mass hypothesis of the 'default' fit of the two particles is the same,
231 //
232 // avoid re-doing the calculation if a pair with the flipped mdst indexes in the map already exists.
233 if (m_useHighestProbMassForExt || (iParticle->getPDGCodeUsedForFit() == jParticle->getPDGCodeUsedForFit())) {
234 if (particleMdstSourcePairsToDist.count({jMdstSource, iMdstSource})) {
235 particleMdstSourcePairsToDist[partMdstSourcePair] = particleMdstSourcePairsToDist[ {jMdstSource, iMdstSource}];
236 continue;
237 }
238 }
239 // Calculate the pair-wise distance.
240 particleMdstSourcePairsToDist[partMdstSourcePair] = this->getDistAtDetSurface(iParticle, jParticle, iDetLayer);
241 }
242
243 }
244
245 // For each particle in the input list, find the minimum among all distances to the reference particles.
246 for (const auto& targetParticle : targetParticles) {
247
248 auto iMdstSource = targetParticle.first;
249 auto iParticle = targetParticle.second;
250
251 // Save the distances and the mdst indexes of the reference particles.
252 std::vector<std::pair<double, int>> iDistancesAndRefMdstSources;
253 for (const auto& [mdstIdxs, dist] : particleMdstSourcePairsToDist) {
254 if (mdstIdxs.first == iMdstSource) {
255 if (!std::isnan(dist) && dist >= 0) {
256 iDistancesAndRefMdstSources.push_back(std::make_pair(dist, mdstIdxs.second));
257 }
258 }
259 }
260
261 if (!iDistancesAndRefMdstSources.size()) {
262 B2DEBUG(12, "The container of distances is empty. Perhaps the target and reference lists contain the same exact particles?");
263 continue;
264 }
265
266 const auto minDist = *std::min_element(iDistancesAndRefMdstSources.begin(), iDistancesAndRefMdstSources.end(),
267 [](const auto & l, const auto & r) {return l.first < r.first;});
268
269 auto jParticle = m_pListReference->getParticleWithMdstSource(minDist.second);
270
271 B2DEBUG(11, "\n"
272 << "Particle w/ mdstSource[" << iMdstSource << "] (PDG = "
273 << iParticle->getPDGCode() << "). Closest charged particle w/ mdstSource["
274 << minDist.second
275 << "] (PDG = " << jParticle->getPDGCode()
276 << ") at " << iDetLayer
277 << " surface is found at D = " << minDist.first
278 << " [cm]\n"
279 << "Storing extraInfo variables:\n"
280 << m_detLayerToDistVariable[iDetLayer]
281 << "\n"
282 << m_detLayerToRefPartIdxVariable[iDetLayer]);
283
284 if (!iParticle->hasExtraInfo(m_detLayerToDistVariable[iDetLayer])) {
285 m_particles[iParticle->getArrayIndex()]->addExtraInfo(m_detLayerToDistVariable[iDetLayer], minDist.first);
286 }
287 m_particles[iParticle->getArrayIndex()]->writeExtraInfo(m_detLayerToRefPartIdxVariable[iDetLayer], minDist.second);
288
289 } // end loop over input particle list.
290
291 } // end loop over detector layers.
292
293 } // end loop over input detectors.
294
295 // Now calculate the isolation score for each target particle.
296 for (const auto& targetParticle : targetParticles) {
297
298 auto iMdstSource = targetParticle.first;
299 auto iParticle = targetParticle.second;
300
301 // Initialise isolation score.
302 double isoScore(0.0);
303 double isoScoreAsWeightedAvg(0.0);
304 double sumWeights(0.0);
305
306 B2DEBUG(11, "Particle w/ mdstSource[" << iMdstSource << "] (PDG = " << iParticle->getPDGCode() << ").");
307
308 for (const auto& detName : m_detNames) {
309
310 auto detWeight = this->getDetectorWeight(iParticle, detName);
311
312 auto weightedSumNonIsoLayers = this->getWeightedSumNonIsoLayers(iParticle, detName, detWeight);
313 auto weightedSumInvDists = this->getWeightedSumInvDists(iParticle, detName, detWeight);
314
315 isoScore += weightedSumNonIsoLayers;
316 isoScoreAsWeightedAvg += weightedSumInvDists;
317 sumWeights += detWeight;
318
319 }
320
321 // Normalise the isolation score to lie in [0, 1].
322 auto minScore = 0.;
323 auto maxScore = m_detNames.size();
324 isoScore = (isoScore - minScore) / (maxScore - minScore);
325
326 isoScoreAsWeightedAvg /= sumWeights;
327 // Normalise weighted average between [0, 1].
328 // But first, clip values that are too large.
329 // Finally, ensure a larger score is assigned for well-isolated tracks.
330 isoScoreAsWeightedAvg = 1. - ((std::min(isoScoreAsWeightedAvg, 1e2) - minScore) / (1e2 - minScore));
331
332 B2DEBUG(11, "\n"
333 << "Isolation score: " << isoScore
334 << "\n"
335 << "Isolation score (as weighted avg): " << isoScoreAsWeightedAvg
336 << "\n"
337 << "Storing extraInfo variable:\n"
339 << "\n"
341
342 m_particles[iParticle->getArrayIndex()]->writeExtraInfo(m_isoScoreVariable, isoScore);
343 m_particles[iParticle->getArrayIndex()]->writeExtraInfo(m_isoScoreVariableAsWeightedAvg, isoScoreAsWeightedAvg);
344
345 }
346
347 B2DEBUG(11, "Finished processing EVENT: " << m_event_metadata->getEvent());
348
349}
350
351
352double TrackIsoCalculatorModule::getDetectorWeight(const Particle* iParticle, const std::string& detName) const
353{
354
356 return -1.;
357 }
358
359 auto hypo = Const::ChargedStable(std::abs(iParticle->getPDGCode()));
360 auto p = iParticle->getP();
361 auto theta = iParticle->getMomentum().Theta();
362 auto det = this->getDetEnum(detName);
363
364 auto detWeight = (*m_DBWeights.get())->getWeight(hypo, det, p, theta);
365 // If w > 0, the detector has detrimental impact on PID:
366 // the value is reset to zero to prevent the detector from contributing to the score.
367 // NB: NaN should stay NaN.
368 detWeight = (detWeight < 0 || std::isnan(detWeight)) ? detWeight : 0.0;
369
370 return detWeight;
371
372}
373
374
375double TrackIsoCalculatorModule::getWeightedSumNonIsoLayers(const Particle* iParticle, const std::string& detName,
376 const float detWeight) const
377{
378
379 auto det = this->getDetEnum(detName);
380 auto nLayers = DetectorSurface::detToLayers.at(detName).size();
381
382 unsigned int n(0);
383 for (const auto& iLayer : DetectorSurface::detToLayers.at(detName)) {
384 auto iDetLayer = detName + std::to_string(iLayer);
385 auto distVar = m_detLayerToDistVariable.at(iDetLayer);
386 auto threshDist = this->getDistThreshold(det, iLayer);
387
388 if (iParticle->hasExtraInfo(distVar)) {
389 if (iParticle->getExtraInfo(distVar) < threshDist) {
390 n++;
391 }
392 }
393
394 }
395
396 if (!n) {
397 B2DEBUG(12, "\nNo close-enough neighbours to this particle in the " << detName << " were found.");
398 }
399
400 if (n > nLayers) {
401 B2FATAL("\nTotal layers in " << detName << ": N=" << nLayers
402 << "\n"
403 << "Layers w/ a close-enough particle: n=" << n
404 << "\n"
405 << "n > N ?? Abort...");
406 }
407
408 return 1. - (-detWeight * (n / nLayers));
409
410}
411
412
413double TrackIsoCalculatorModule::getWeightedSumInvDists(const Particle* iParticle, const std::string& detName,
414 const float detWeight) const
415{
416
417 auto det = this->getDetEnum(detName);
418
419 double sumInvDists(0.);
420 for (const auto& iLayer : DetectorSurface::detToLayers.at(detName)) {
421
422 auto iDetLayer = detName + std::to_string(iLayer);
423 auto distVar = m_detLayerToDistVariable.at(iDetLayer);
424 auto threshDist = this->getDistThreshold(det, iLayer);
425
426 if (iParticle->hasExtraInfo(distVar)) {
427 sumInvDists += threshDist / iParticle->getExtraInfo(distVar);
428 }
429
430 }
431
432 return detWeight * sumInvDists;
433
434}
435
436
438 const Particle* jParticle,
439 const std::string& detLayerName) const
440{
441
442 // Radius and z boundaries of the cylinder describing this detector's surface.
443 const auto rho = DetectorSurface::detLayerToSurfBoundaries.at(detLayerName).m_rho;
444 const auto zfwd = DetectorSurface::detLayerToSurfBoundaries.at(detLayerName).m_zfwd;
445 const auto zbwd = DetectorSurface::detLayerToSurfBoundaries.at(detLayerName).m_zbwd;
446 // Polar angle boundaries between barrel and endcaps.
447 const auto th_fwd = DetectorSurface::detLayerToSurfBoundaries.at(detLayerName).m_th_fwd;
448 const auto th_fwd_brl = DetectorSurface::detLayerToSurfBoundaries.at(detLayerName).m_th_fwd_brl;
449 const auto th_bwd_brl = DetectorSurface::detLayerToSurfBoundaries.at(detLayerName).m_th_bwd_brl;
450 const auto th_bwd = DetectorSurface::detLayerToSurfBoundaries.at(detLayerName).m_th_bwd;
451
452 std::string nameExtTheta = "helixExtTheta(" + std::to_string(rho) + "," + std::to_string(zfwd) + "," + std::to_string(zbwd) + ")";
454 nameExtTheta.replace(nameExtTheta.size() - 1, 1, ", 1)");
455 }
456 const auto iExtTheta = std::get<double>(Variable::Manager::Instance().getVariable(nameExtTheta)->function(iParticle));
457 const auto jExtTheta = std::get<double>(Variable::Manager::Instance().getVariable(nameExtTheta)->function(jParticle));
458
459 std::string nameExtPhi = "helixExtPhi(" + std::to_string(rho) + "," + std::to_string(zfwd) + "," + std::to_string(zbwd) + ")";
461 nameExtPhi.replace(nameExtPhi.size() - 1, 1, ", 1)");
462 }
463 const auto iExtPhi = std::get<double>(Variable::Manager::Instance().getVariable(nameExtPhi)->function(iParticle));
464 const auto jExtPhi = std::get<double>(Variable::Manager::Instance().getVariable(nameExtPhi)->function(jParticle));
465
466 const auto iExtInBarrel = (iExtTheta >= th_fwd_brl && iExtTheta < th_bwd_brl);
467 const auto jExtInBarrel = (jExtTheta >= th_fwd_brl && jExtTheta < th_bwd_brl);
468
469 const auto iExtInFWDEndcap = (iExtTheta >= th_fwd && iExtTheta < th_fwd_brl);
470 const auto jExtInFWDEndcap = (jExtTheta >= th_fwd && jExtTheta < th_fwd_brl);
471
472 const auto iExtInBWDEndcap = (iExtTheta >= th_bwd_brl && iExtTheta < th_bwd);
473 const auto jExtInBWDEndcap = (jExtTheta >= th_bwd_brl && jExtTheta < th_bwd);
474
475 if (boost::contains(detLayerName, "CDC") || boost::contains(detLayerName, "TOP")) {
476
477 // If any of the two extrapolated tracks is not in the barrel region of the CDC/TOP, the distance is undefined.
478 if (!iExtInBarrel || !jExtInBarrel) {
479 return std::numeric_limits<double>::quiet_NaN();
480 }
481
482 // For CDC and TOP, we calculate the distance between the points where the two input
483 // extrapolated track helices cross the input detector's cylindrical surface
484 // on the (rho, phi) plane. Namely, this is the cord length of the arc
485 // that subtends deltaPhi.
486 auto diffPhi = jExtPhi - iExtPhi;
487 if (std::abs(diffPhi) > M_PI) {
488 diffPhi = (diffPhi > M_PI) ? diffPhi - 2 * M_PI : 2 * M_PI + diffPhi;
489 }
490
491 return 2.0 * rho * sin(std::abs(diffPhi) / 2.0);
492
493 } else {
494
495 if (boost::contains(detLayerName, "ARICH")) {
496 // If any of the two tracks is not in the ARICH theta acceptance, the distance is undefined.
497 if (!iExtInFWDEndcap || !jExtInFWDEndcap) {
498 return std::numeric_limits<double>::quiet_NaN();
499 }
500 }
501 if (boost::contains(detLayerName, "ECL") || boost::contains(detLayerName, "KLM")) {
502
503 // For ECL and KLM, we require track pairs to be both in the barrel,
504 // both in the FWD endcap, or both in the FWD endcap. Otherwise, the distance is undefined.
505 if (
506 !(iExtInBarrel && jExtInBarrel) &&
507 !(iExtInFWDEndcap && jExtInFWDEndcap) &&
508 !(iExtInBWDEndcap && jExtInBWDEndcap)
509 ) {
510 return std::numeric_limits<double>::quiet_NaN();
511 }
512 }
513
514 // Ok, we know theta and phi.
515 // Let's extract (spherical) R by using the transformation:
516 //
517 // 1. rho = r * sin(theta)
518 // 2. phi = phi
519 // 3. z = r * cos(theta)
520 //
521 // The formula to be inverted depends on where each extrapolated track's theta is found:
522 // if in barrel, use 1. (rho is known), if in fwd/bwd endcap, use 3. (zfwd/zbwd is known).
523
524 const auto iExtR = (iExtTheta >= th_fwd_brl && iExtTheta < th_bwd_brl) ? rho / sin(iExtTheta) : ((iExtTheta >= th_fwd
525 && iExtTheta < th_fwd_brl) ? zfwd / cos(iExtTheta) : zbwd / cos(iExtTheta));
526 const auto jExtR = (jExtTheta >= th_fwd_brl && jExtTheta < th_bwd_brl) ? rho / sin(jExtTheta) : ((jExtTheta >= th_fwd
527 && jExtTheta < th_fwd_brl) ? zfwd / cos(jExtTheta) : zbwd / cos(jExtTheta));
528
529 return sqrt((iExtR * iExtR) + (jExtR * jExtR) - 2 * iExtR * jExtR * (sin(iExtTheta) * sin(jExtTheta) * cos(iExtPhi - jExtPhi)
530 + cos(iExtTheta) * cos(jExtTheta)));
531
532 }
533
534 return std::numeric_limits<double>::quiet_NaN();
535
536}
537
539{
540
541 bool checkPList = false;
542
544 checkPList = Const::chargedStableSet.find(abs(m_decaydescriptor.getMother()->getPDGCode())) != Const::invalidParticle;
545 } else {
546 for (int pdgcode : m_decaydescriptor.getSelectionPDGCodes()) {
547 checkPList = Const::chargedStableSet.find(abs(pdgcode)) != Const::invalidParticle;
548 if (!checkPList) {
549 break;
550 }
551 }
552 }
553
555 bool checkPListReference = false;
556 if (dd.init(m_pListReferenceName)) {
557 checkPListReference = Const::chargedStableSet.find(abs(dd.getMother()->getPDGCode())) != Const::invalidParticle;
558 }
559
560 return (checkPList and checkPListReference);
561
562};
Provides a type-safe way to pass members of the chargedStableSet set.
Definition Const.h:589
static const ParticleSet chargedStableSet
set of charged stable particles
Definition Const.h:618
static const ParticleType invalidParticle
Invalid particle, used internally.
Definition Const.h:681
Represents a particle in the DecayDescriptor.
int getPDGCode() const
Return PDG code.
std::string getFullName() const
returns the full name of the particle full_name = name:label
The DecayDescriptor stores information about a decay tree or parts of a decay tree.
bool init(const std::string &str)
Initialise the DecayDescriptor from given string.
const DecayDescriptorParticle * getMother() const
return mother.
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition Module.cc:208
Module()
Constructor.
Definition Module.cc:30
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition Module.h:80
Class to store reconstructed particles.
Definition Particle.h:76
bool hasExtraInfo(const std::string &name) const
Return whether the extra info with the given name is set.
Definition Particle.cc:1351
int getPDGCode(void) const
Returns PDG code.
Definition Particle.h:465
ROOT::Math::XYZVector getMomentum() const
Returns momentum vector.
Definition Particle.h:580
double getP() const
Returns momentum magnitude (same as getMomentumMagnitude but with shorter name)
Definition Particle.h:598
double getExtraInfo(const std::string &name) const
Return given value if set.
Definition Particle.cc:1374
StoreObjPtr< EventMetaData > m_event_metadata
The event information.
std::unique_ptr< DBObjPtr< PIDDetectorWeights > > m_DBWeights
Interface to get the database payload with the PID detector weights.
~TrackIsoCalculatorModule() override
Destructor, use this to clean up anything you created in the constructor.
void initialize() override
Use this to initialize resources or memory your module needs.
double getDistAtDetSurface(const Particle *iParticle, const Particle *jParticle, const std::string &detLayerName) const
Calculate the distance between the points where the two input extrapolated track helices cross the gi...
void event() override
Called once for each event.
std::unordered_map< std::string, std::string > m_detLayerToDistVariable
Map that associates to each detector layer (e.g., 'CDC6') the name of the variable representing the d...
std::string m_decayString
The name of the input charged stable particle list, or composite particle w/ charged stable daughters...
StoreObjPtr< ParticleList > m_pListTarget
The input ParticleList object for which distances are to be calculated.
std::string m_isoScoreVariable
The name of the variable representing the track isolation score.
double getWeightedSumInvDists(const Particle *iParticle, const std::string &detName, const float detWeight) const
Get the sum of the inverse (scaled) minimum distances over the given detector layers,...
StoreArray< Particle > m_particles
StoreArray of Particles.
std::string m_payloadName
The name of the database payload object with the MVA weights.
std::vector< std::string > m_detNames
The list of names of the detectors at whose inner (cylindrical) surface we extrapolate each track's p...
Const::EDetector getDetEnum(const std::string &detName) const
Get the enum type for this detector name.
std::string m_pListReferenceName
The name of the input ParticleList of reference tracks.
std::unordered_map< std::string, std::string > m_detLayerToRefPartIdxVariable
Map that associates to each detector layer (e.g, 'CDC6') the name of the variable representing the md...
double getDistThreshold(Const::EDetector det, int layer) const
Get the threshold value per detector layer for the distance to closest ext.
bool onlySelectedStdChargedInDecay()
Check whether input particle list and reference list are of a valid charged stable particle.
double getDetectorWeight(const Particle *iParticle, const std::string &detName) const
Get the PID weight, , for this particle and detector reading it from the payload, if selected.
bool m_excludePIDDetWeights
Exclude the PID detector weights for the isolation score definition.
DecayDescriptor m_decaydescriptor
< Decay descriptor of decays to look for.
TrackIsoCalculatorModule()
Constructor: Sets the description, the properties and the parameters of the module.
double getWeightedSumNonIsoLayers(const Particle *iParticle, const std::string &detName, const float detWeight) const
Get the sum of layers with a close-by track, divided by the total number of layers,...
unsigned short m_nSelectedDaughters
The number of selected daughters in the decay string.
StoreObjPtr< ParticleList > m_pListReference
The input ParticleList object of reference tracks.
std::string m_isoScoreVariableAsWeightedAvg
The name of the variable representing the track isolation score.
bool m_useHighestProbMassForExt
If this option is set, the helix extrapolation for the target and reference particles will use the tr...
static Manager & Instance()
get singleton instance.
Definition Manager.cc:26
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition Module.h:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
double sqrt(double a)
sqrt for double
Definition beamHelpers.h:28
Abstract base class for different kinds of events.
static const std::unordered_map< std::string, DetSurfCylBoundaries > detLayerToSurfBoundaries
Map that associates to each detector layer its valid cylindrical surface's boundaries.
static const std::unordered_map< std::string, std::vector< int > > detToLayers
Map that associates to each detector its list of valid layers.