PYTHIA  8.312
ResonanceDecayFilterHook.h
1 // ResonanceDecayFilterHook.h is part of the PYTHIA event generator.
2 // Copyright (C) 2024 Stephen Mrenna, Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL v2 or later, see COPYING for details.
4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
5 // Author: Stephen Mrenna, January 2021.
6 
7 // This class is used to filter out resonance decay products to obtain
8 // a desired parton-level final state. It copies an implementation
9 // already in the CMS software originally authored by Josh Bendavid,
10 // with some minor modifications.
11 
12 // The following settings are available:
13 // flag ResonanceDecayFilter:filter = false
14 // Activate filter if true.
15 // flag ResonanceDecayFilter:exclusive = false
16 // Demand exactly the particles requested by the filter if true.
17 // If false, additional particles can be present in the resonance decays.
18 // flag ResonanceDecayFilter:eMuAsEquivalent = false
19 // If true, treat the particle id for electron and muon as the same.
20 // flag ResonanceDecayFilter:eMuTauAsEquivalent = false
21 // As above, but include tau leptons.
22 // flag ResonanceDecayFilter:allNuAsEquivalent = false
23 // flag ResonanceDecayFilter:udscAsEquivalent = false
24 // If true, treat the particle id for all up, down, strange, and charm
25 // quarks as the same.</flag>
26 // flag ResonanceDecayFilter:udscbAsEquivalent = false
27 // As above, but include bottom quarks.
28 // flag ResonanceDecayFilter:wzAsEquivalent = false
29 // If true, treat the particle id for W and Z bosons as the same.
30 // mvec ResonanceDecayFilter:mothers = {}
31 // If provided, the set of mothers that will be filtered.
32 // Otherwise, all resonance decay products are considered.
33 // mvec ResonanceDecayFilter:daughters = {}
34 // A list of particle ids for resonance decay products that
35 // will pass the filter. If empty, all particles will pass the
36 // filter. Beware: an infinite loop is possible if the daughter never
37 // appears in resonance decays.
38 
39 #ifndef Pythia8_ResonanceDecayFilterHooks_H
40 #define Pythia8_ResonanceDecayFilterHooks_H
41 
42 // Includes.
43 #include "Pythia8/Pythia.h"
44 #include "Pythia8/UserHooks.h"
45 #include "Pythia8/Event.h"
46 
47 namespace Pythia8 {
48 
50 public:
51 
52  // Constructor.
54 
55  // Override base class methods.
56  bool canVetoResonanceDecays() override {return true;}
57  bool doVetoResonanceDecays(Event& process) override {
58  return checkVetoResonanceDecays(process);}
59  bool initAfterBeams() override;
60 
61  // Class specific.
62  bool checkVetoResonanceDecays(const Event& process);
63  unsigned long int returnCounter() {return counter;};
64 
65 private:
66 
67  // Return the particle ID category.
68  int idCat(int id);
69 
70  // Data members.
71  bool filter, exclusive, eMuAsEquivalent, eMuTauAsEquivalent,
72  allNuAsEquivalent, udscAsEquivalent, udscbAsEquivalent, wzAsEquivalent;
73  unsigned long int counter;
74  set<int> mothers;
75  vector<int> daughters;
76  unordered_map<int, int> requestedDaughters, observedDaughters;
77 
78 };
79 
80 //--------------------------------------------------------------------------
81 
82 // Constructor.
83 
85  counter = 0;
86  settings.addFlag("ResonanceDecayFilter:filter", false);
87  settings.addFlag("ResonanceDecayFilter:exclusive", false);
88  settings.addFlag("ResonanceDecayFilter:eMuAsEquivalent", false);
89  settings.addFlag("ResonanceDecayFilter:eMuTauAsEquivalent", false);
90  settings.addFlag("ResonanceDecayFilter:allNuAsEquivalent", false);
91  settings.addFlag("ResonanceDecayFilter:udscAsEquivalent", false);
92  settings.addFlag("ResonanceDecayFilter:udscbAsEquivalent", false);
93  settings.addFlag("ResonanceDecayFilter:wzAsEquivalent", false);
94  settings.addMVec("ResonanceDecayFilter:mothers", vector<int>(), false,
95  false, 0, 0);
96  settings.addMVec("ResonanceDecayFilter:daughters", vector<int>(), false,
97  false, 0, 0);
98 }
99 
100 //--------------------------------------------------------------------------
101 
102 // Return a particle ID given equivalence user settings.
103 
104 int ResonanceDecayFilterHook::idCat(int id) {
105  id = abs(id);
106  if (id == 13 && (eMuAsEquivalent || eMuTauAsEquivalent)) id = 11;
107  else if (id == 15 && eMuTauAsEquivalent) id = 11;
108  else if ((id == 14 || id == 16) && allNuAsEquivalent) id = 12;
109  else if ((id == 2 || id == 3 || id == 4) && udscAsEquivalent) id = 1;
110  else if ((id == 2 || id == 3 || id == 4 || id == 5) &&
111  udscbAsEquivalent) id = 1;
112  else if ((id == 23 || id == 24) && wzAsEquivalent) id = 23;
113  return id;
114 }
115 
116 //--------------------------------------------------------------------------
117 
118 // Intialize the user hook after the beams.
119 
121  filter = settingsPtr->flag("ResonanceDecayFilter:filter");
122  exclusive = settingsPtr->flag("ResonanceDecayFilter:exclusive");
123  eMuAsEquivalent = settingsPtr->flag("ResonanceDecayFilter:eMuAsEquivalent");
124  eMuTauAsEquivalent = settingsPtr->
125  flag("ResonanceDecayFilter:eMuTauAsEquivalent");
126  allNuAsEquivalent = settingsPtr->flag
127  ("ResonanceDecayFilter:allNuAsEquivalent");
128  udscAsEquivalent = settingsPtr->
129  flag("ResonanceDecayFilter:udscAsEquivalent");
130  udscbAsEquivalent = settingsPtr->
131  flag("ResonanceDecayFilter:udscbAsEquivalent");
132  wzAsEquivalent = settingsPtr->flag("ResonanceDecayFilter:wzAsEquivalent");
133  auto mothersIn = settingsPtr->mvec("ResonanceDecayFilter:mothers");
134  mothers.clear();
135  mothers.insert(mothersIn.begin(), mothersIn.end());
136  daughters = settingsPtr->mvec("ResonanceDecayFilter:daughters");
137  requestedDaughters.clear();
138 
139  // Loop over the daughters.
140  for (int id : daughters) ++requestedDaughters[idCat(id)];
141  return true;
142 
143 }
144 
145 //--------------------------------------------------------------------------
146 
147 // Return true of the resonance decays are vetoed.
148 
150  if (!filter) return false;
151 
152  // Count the number of times hook is called.
153  counter++;
154  observedDaughters.clear();
155 
156  // Loop over particles and determine equivalent types.
157  for (int i = 0; i < process.size(); ++i) {
158  const Particle &p = process[i];
159  int mid = p.mother1() > 0 ? abs(process[p.mother1()].id()) : 0;
160 
161  // If no list of mothers is provided, then all particles
162  // in hard process and resonance decays are counted together
163  if (mothers.empty() || mothers.count(mid) || mothers.count(-mid))
164  ++observedDaughters[idCat(p.id())];
165  }
166 
167  // Check if criteria is satisfied.
168  // inclusive mode: at least as many decay products as requested
169  // exclusive mode: exactly as many decay products as requested
170  // (but additional particle types not appearing in the list of requested
171  // daughter id's are ignored)
172  for (const auto &reqDau : requestedDaughters) {
173  int reqId = reqDau.first;
174  int reqCount = reqDau.second;
175  auto obsItr = observedDaughters.find(reqId);
176  int obsCount = ( obsItr != observedDaughters.end() ) ? obsItr->second : 0;
177 
178  // Inclusive criteria not satisfied, veto event
179  if (obsCount < reqCount) return true;
180 
181  // Exclusive criteria not satisfied, veto event
182  if (exclusive && obsCount > reqCount) return true;
183  }
184 
185  // All criteria satisfied, don't veto
186  return false;
187 
188 }
189 
190 //==========================================================================
191 
192 } // end namespace Pythia8
193 
194 #endif // end Pythia8_ResonanceDecayHooks_H
bool flag(string keyIn)
Give back current value, with check that key exists.
Definition: Settings.cc:1478
Settings * settingsPtr
Pointer to the settings database.
Definition: PhysicsBase.h:81
The Event class holds all info on the generated event.
Definition: Event.h:453
ResonanceDecayFilterHook(Settings &settings)
Constructor.
Definition: ResonanceDecayFilterHook.h:84
void id(int idIn)
Member functions for input.
Definition: Event.h:88
bool initAfterBeams() override
Intialize the user hook after the beams.
Definition: ResonanceDecayFilterHook.h:120
Definition: ResonanceDecayFilterHook.h:49
int size() const
Event record size.
Definition: Event.h:504
Definition: Event.h:32
bool canVetoResonanceDecays() override
Override base class methods.
Definition: ResonanceDecayFilterHook.h:56
bool doVetoResonanceDecays(Event &process) override
Definition: ResonanceDecayFilterHook.h:57
Header for classes to set beam momentum and interaction vertex spread.
Definition: Analysis.h:20
UserHooks is base class for user access to program execution.
Definition: UserHooks.h:32
bool checkVetoResonanceDecays(const Event &process)
Class specific.
Definition: ResonanceDecayFilterHook.h:149
bool flag(string key) const
Shorthand to read settings values.
Definition: PhysicsBase.h:45
Definition: Settings.h:195
void addFlag(string keyIn, bool defaultIn)
Add new entry.
Definition: Settings.h:267