ATLAS Offline Software
DecayToSUEP.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2020-2022 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 #include "UserHooksUtils.h"
6 #include "UserSetting.h"
8 #include <boost/phoenix/bind.hpp>
9 #include <boost/math/tools/roots.hpp>
10 #include <math.h>
11 #include <stdexcept>
12 #include <iostream>
13 
14 #include "suep_shower.h"
15 
16 namespace Pythia8{
17  class DecayToSUEP;
18 }
19 
21 
22 // Uncomment the following line to enable debug messages to be printed to std::cout
23 //#define SUEP_DEBUG 1
24 
25 namespace Pythia8{
26 
27  //******************************
28  //** Main UserHook derived class
29  //******************************
30 
39  class DecayToSUEP : public UserHooks{
40 
41  public:
42 
44  m_pdgId("DecayToSUEP:PDGId", 25),
45  m_mass("DecayToSUEP:Mass", 125.0),
46  m_darkMesonMass("DecayToSUEP:DarkMesonMass", 1.),
47  m_darkTemperature("DecayToSUEP:DarkTemperature", 1.) {
48 
49  std::cout<<"**********************************************************"<<std::endl;
50  std::cout<<"* *"<<std::endl;
51  std::cout<<"* Enabled SUEP decay UserHook! *"<<std::endl;
52  std::cout<<"* *"<<std::endl;
53  std::cout<<"**********************************************************"<<std::endl;
54 
55  }
56 
58 
59  // Enable the call to the user-hook
60  virtual bool canVetoProcessLevel() {
61  return true;
62  }
63 
64  /* Actually implement the user hook.
65  *
66  * We slightly abuse this function in the sense that no veto is performed but instead
67  * the hook is used to modify the event record in between the process-level and parton-level steps.
68  * This modification is allowed in the Pythia8 manual, with warning of event consistency that
69  * should pose no harm in this context.
70  */
71  virtual bool doVetoProcessLevel(Event& process);
72 
73 
74  private:
75 
78 
81 
84 
87 
88  };
89 
90 
91  //******************************
92  //** Implementation of DecayToSUEP UserHook
93  //******************************
95  {
96 #ifdef SUEP_DEBUG
97  std::cout << "[SUEP_DEBUG] " << "Start of user hook for this event." << std::endl;
98 #endif
99 
100  //First, find the particle to decay
101  bool particleFound=false;
102 
103 #ifdef SUEP_DEBUG
104  for (int ii=0; ii < process.size(); ++ii) {
105  std::cout << "[SUEP_DEBUG] " << ii << ": id=" << process[ii].id() << ", Final=" << process[ii].isFinal() << ", Status=" << process[ii].status() << ", daughter1=" << process[ii].daughter1() << ", daughter2=" << process[ii].daughter2() << std::endl;
106  }
107 #endif
108 
109  for (int ii=0; ii < process.size(); ++ii) {
110  if ( (process[ii].id() == m_pdgId(settingsPtr)) and (process[ii].daughter1()!=process[ii].daughter2() && process[ii].daughter1()>0 && process[ii].daughter2()>0) ) {
111 
112  Vec4 higgs4mom, mesonmom;
113  vector< Vec4 > suep_shower4momenta;
114  particleFound=true;
115 
116  //setup SUEP shower
117  static Suep_shower suep_shower(m_darkMesonMass(settingsPtr), m_darkTemperature(settingsPtr), m_mass(settingsPtr), rndmPtr);
118 
119 #ifdef SUEP_DEBUG
120  std::cout << "[SUEP_DEBUG] " << "Particle (pdgId=" << m_pdgId(settingsPtr) << ", isFinal=True) found. Decaying to SUEP now." << std::endl;
121 #endif
122 
123  // First undo decay
124  process[ii].undoDecay();
125 
126  int originalEventSize = process.size();
127 
128  // Generate the shower, output are 4 vectors in the rest frame of the shower
129  higgs4mom=process[ii].p();
130  int nShowerAttempts=0;
131  do {
132  try {
133  nShowerAttempts++;
134  suep_shower4momenta=suep_shower.generate_shower();
135  if( suep_shower4momenta.size()<3){
136  //Failed to balance energy or less than 3 particles in the shower
137  //Try again until nShowerAttempts >= 3
138  } else {
139  //All ok!
140  nShowerAttempts = -1; //exit condition
141  }
142  } catch (std::exception &e) {
143  //Failed to generate the shower!
144  //Can happen in some rare circumstances,
145  //Try again until nShowerAttempts >= 3
146  }
147  } while ((nShowerAttempts > 0) && (nShowerAttempts < 3));
148  if (nShowerAttempts >= 3) {
149  //Something is seriously wrong then, print warning and skip to next event
150  std::cout << "[SUEP] WARNING: Something went wrong in generating the shower. Skipping the event." << std::endl;
151  return true; //veto the event!
152  }
153 
154  // Loop over hidden sector mesons and append to the event
155  for (unsigned j = 0; j < suep_shower4momenta.size(); ++j){
156  //construct pythia 4vector
157  mesonmom = suep_shower4momenta[j];
158 
159  // boost to the lab frame
160  mesonmom.bst(higgs4mom.px()/higgs4mom.e(),higgs4mom.py()/higgs4mom.e(), higgs4mom.pz()/higgs4mom.e());
161 
162  //append particle to the event. Hidden/dark meson pdg code is 999999.
163  process.append(999999, 91, ii, 0, 0, 0, 0, 0, mesonmom.px(), mesonmom.py(), mesonmom.pz(), mesonmom.e(), m_darkMesonMass(settingsPtr));
164 
165 #ifdef SUEP_DEBUG
166  std::cout << "[SUEP_DEBUG] " << "Adding dark meson with px=" << mesonmom.px() << ", py=" << mesonmom.py() << ", pz=" << mesonmom.pz() << ", m=" << m_darkMesonMass(settingsPtr) << std::endl;
167 #endif
168  }
169 
170  // Just to be sure, only modify Higgs status and daughters if a valid decay did happen
171  if ( suep_shower4momenta.size() > 0 ) {
172 #ifdef SUEP_DEBUG
173  std::cout << "[SUEP_DEBUG] " << "Setting original particle status-code as non-Final particle. Adding daughters with indices: " << originalEventSize << " - " << process.size()-1 << std::endl;
174 #endif
175  // Change the status code of the Higgs to reflect that it has decayed.
176  process[ii].statusNeg();
177 
178  //set daughters of the Higgs. Take advantage that we just appended them
179  process[ii].daughters(originalEventSize, process.size()-1);
180  }
181 
182  //no need to continue the loop
183  break;
184 
185  } // if particle to decay found
186 
187  } // loop over particles in the event
188 
189  if (not particleFound) {
190  std::cout << "[DecayToSUEP] " << "Particle " << m_pdgId(settingsPtr) << " not found. Nothing to decay to SUEP for this event." << std::endl;
191  } else {
192 #ifdef SUEP_DEBUG
193  std::cout << "[SEUP_DEBUG] " << "All Done for this event." << std::endl;
194 #endif
195  }
196 #ifdef SUEP_DEBUG
197  std::cout << "[SUEP_DEBUG] Printing event after adding SUEP:" << std::endl;
198  for (int ii=0; ii < process.size(); ++ii) {
199  std::cout << "[SUEP_DEBUG] " << ii << ": id=" << process[ii].id() << ", Final=" << process[ii].isFinal() << ", mayDecay=" << process[ii].mayDecay() << ", Status=" << process[ii].status() << ", daughter1=" << process[ii].daughter1() << ", daughter2=" << process[ii].daughter2() << std::endl;
200  }
201 #endif
202 
203  //return false: let the event continue
204  return false;
205  }
206 
207 
208 } // namespace Pythia8
Pythia8::DecayToSUEP
Pythia8 UserHook to decay a given scalar particle to SUEP.
Definition: DecayToSUEP.cxx:39
AllowedVariables::e
e
Definition: AsgElectronSelectorTool.cxx:37
Suep_shower::generate_shower
std::vector< Pythia8::Vec4 > generate_shower()
Generate a shower event, in the rest frame of the showe.
Definition: suep_shower.cxx:123
Event
Definition: trigbs_orderedMerge.cxx:42
DecayToSUEPCreator
Pythia8_UserHooks::UserHooksFactory::Creator< Pythia8::DecayToSUEP > DecayToSUEPCreator("DecayToSUEP")
Pythia8_UserHooks::UserHooksFactory::Creator
Definition: UserHooksFactory.h:54
Pythia8::DecayToSUEP::m_darkTemperature
Pythia8_UserHooks::UserSetting< double > m_darkTemperature
Temperature parameter [GeV].
Definition: DecayToSUEP.cxx:86
SUSY_SimplifiedModel_PostInclude.process
string process
Definition: SUSY_SimplifiedModel_PostInclude.py:42
Pythia8::DecayToSUEP::~DecayToSUEP
~DecayToSUEP()
Definition: DecayToSUEP.cxx:57
UserHooksFactory.h
Pythia8
Author: James Monk (jmonk@cern.ch)
Definition: IPythia8Custom.h:13
Pythia8_UserHooks::UserSetting< int >
suep_shower.h
Suep_shower
Auxiliary class for SUEP generation.
Definition: suep_shower.h:33
Pythia8::DecayToSUEP::m_mass
Pythia8_UserHooks::UserSetting< double > m_mass
Mass of system decaying [GeV].
Definition: DecayToSUEP.cxx:80
Pythia8::DecayToSUEP::DecayToSUEP
DecayToSUEP()
Definition: DecayToSUEP.cxx:43
calibdata.exception
exception
Definition: calibdata.py:496
Pythia8::DecayToSUEP::m_pdgId
Pythia8_UserHooks::UserSetting< int > m_pdgId
PDG Id of particle to be decayed to SUEP.
Definition: DecayToSUEP.cxx:77
UserSetting.h
Pythia8::DecayToSUEP::canVetoProcessLevel
virtual bool canVetoProcessLevel()
Definition: DecayToSUEP.cxx:60
UserHooksUtils.h
Pythia8::DecayToSUEP::m_darkMesonMass
Pythia8_UserHooks::UserSetting< double > m_darkMesonMass
Dark-meson mass parameter [GeV].
Definition: DecayToSUEP.cxx:83
Pythia8::DecayToSUEP::doVetoProcessLevel
virtual bool doVetoProcessLevel(Event &process)
Definition: DecayToSUEP.cxx:94