PYTHIA  8.312
Dire.h
1 // Dire.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2024 Stefan Prestel, 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 
6 // Header file for the top level Dire class.
7 
8 #ifndef Pythia8_Dire_H
9 #define Pythia8_Dire_H
10 
11 #define DIRE_VERSION "2.002"
12 
13 // DIRE includes.
14 #include "Pythia8/DireSplittingLibrary.h"
15 #include "Pythia8/DireMerging.h"
16 #include "Pythia8/DireTimes.h"
17 #include "Pythia8/DireSpace.h"
18 #include "Pythia8/DireWeightContainer.h"
19 #include "Pythia8/DireHooks.h"
20 
21 // Pythia includes.
22 #include "Pythia8/Info.h"
23 #include "Pythia8/Settings.h"
24 #include "Pythia8/ParticleData.h"
25 #include "Pythia8/Basics.h"
26 #include "Pythia8/PartonSystems.h"
27 #include "Pythia8/UserHooks.h"
28 #include "Pythia8/MergingHooks.h"
29 #include "Pythia8/PartonVertex.h"
30 #include "Pythia8/StandardModel.h"
31 #include "Pythia8/ShowerModel.h"
32 
33 #include <iostream>
34 #include <sstream>
35 
36 namespace Pythia8 {
37 
38 //==========================================================================
39 
40 class Dire : public ShowerModel {
41 
42  public:
43 
44  Dire() : weightsPtr(nullptr), timesPtr(nullptr), timesDecPtr(nullptr),
45  spacePtr(nullptr), splittings(nullptr), hooksPtr(nullptr),
46  mergingPtr(nullptr), hardProcessPtr(nullptr), mergingHooksPtr(nullptr),
47  hasOwnWeights(false), hasOwnTimes(false), hasOwnTimesDec(false),
48  hasOwnSpace(false), hasOwnSplittings(false), hasOwnHooks(false),
49  hasUserHooks(false), hasOwnHardProcess(false),
50  hasOwnMergingHooks(false), initNewSettings(false), isInit(false),
51  isInitShower(false), printBannerSave(true) { createPointers(); }
52 
53  Dire( MergingHooksPtr mergingHooksPtrIn, PartonVertexPtr partonVertexPtrIn)
54  : pythiaMergingHooksPtr(mergingHooksPtrIn),
55  partonVertexPtr(partonVertexPtrIn), weightsPtr(nullptr),
56  timesPtr(nullptr), timesDecPtr(nullptr), spacePtr(nullptr),
57  splittings(nullptr), hooksPtr(nullptr),
58  mergingPtr(nullptr), hardProcessPtr(nullptr),
59  hasOwnWeights(false), hasOwnTimes(false), hasOwnTimesDec(false),
60  hasOwnSpace(false), hasOwnSplittings(false), hasOwnHooks(false),
61  hasUserHooks(false), hasOwnHardProcess(false),
62  hasOwnMergingHooks(false), initNewSettings(false), isInit(false),
63  isInitShower(false), printBannerSave(true) { createPointers(); }
64 
65  ~Dire() {
66  if (hasOwnWeights) delete weightsPtr;
67  if (hasOwnSplittings) delete splittings;
68  if (hasOwnHardProcess) delete hardProcessPtr;
69  }
70 
71  // Flexible-use call at the beginning of each event in pythia.next().
72  // Currently not used, but should be used for clearing some internal
73  // bookkeeping that is otherewise reset in shower prepare functions.
74  void onBeginEvent() override {
75  return;
76  }
77 
78  // Flexible-use call at the end of each event in pythia.next().
79  // Currently only to accumulate shower weights.
80  void onEndEvent(PhysicsBase::Status status) override {
81  // No finalize in case of failure.
82  if (status == INCOMPLETE) return;
83  // Update the event weight by the Dire shower weight when relevant.
84  // Retrieve the shower weight.
85  weightsPtr->calcWeight(0.);
86  weightsPtr->reset();
87  double pswt = weightsPtr->getShowerWeight();
88  // Multiply the shower weight to the event weight.
89  double wt = infoPtr->weight();
90  infoPtr->weightContainerPtr->setWeightNominal(wt * pswt);
91  }
92 
93  void createPointers();
94 
95  // Initialization function called before beams are set up.
96  // Currently only to register objects as PhysicsBase (=initialize ptrs).
97  bool init(MergingPtr, MergingHooksPtr, PartonVertexPtr, WeightContainer*)
98  override {
99  subObjects.clear();
100  if (mergingHooksPtr) {
101  registerSubObject(*mergingHooksPtr);
102  }
103  if (mergingPtr) {
105  }
106  if (timesPtr) registerSubObject(*timesPtr);
107  if (timesDecPtr) registerSubObject(*timesDecPtr);
108  if (spacePtr) registerSubObject(*spacePtr);
109  return true;
110  }
111 
112  // Initialization function called after beams are set up, used as main
113  // initialization.
114  bool initAfterBeams() override;
115 
116  void initTune();
117  void initShowersAndWeights();
118  void setup(BeamParticle* beamA, BeamParticle* beamB);
119  void printBanner();
120 
121  TimeShowerPtr getTimeShower() const override { return timesPtr; }
122  TimeShowerPtr getTimeDecShower() const override { return timesDecPtr; }
123  SpaceShowerPtr getSpaceShower() const override { return spacePtr; }
124  MergingHooksPtr getMergingHooks() const override { return mergingHooksPtr; }
125  MergingPtr getMerging() const override { return mergingPtr; }
126 
127  MergingHooksPtr pythiaMergingHooksPtr;
128  PartonVertexPtr partonVertexPtr;
129 
130  DireWeightContainer* weightsPtr;
131  shared_ptr<DireTimes> timesPtr;
132  shared_ptr<DireTimes> timesDecPtr;
133  shared_ptr<DireSpace> spacePtr;
134  DireSplittingLibrary* splittings;
135  DireHooks* hooksPtr;
136 
137  DireInfo direInfo;
138 
139  // Pointer to Dire merging objects.
140  shared_ptr<DireMerging> mergingPtr;
141  DireHardProcess* hardProcessPtr;
142  shared_ptr<DireMergingHooks> mergingHooksPtr;
143 
144  bool hasOwnWeights, hasOwnTimes, hasOwnTimesDec, hasOwnSpace,
145  hasOwnSplittings, hasOwnHooks, hasUserHooks,
146  hasOwnHardProcess, hasOwnMergingHooks;
147  bool initNewSettings, isInit, isInitShower, printBannerSave;
148 
149 };
150 
151 //==========================================================================
152 
153 } // end namespace Pythia8
154 
155 #endif // Pythia8_Dire_H
shared_ptr< DireMerging > mergingPtr
Pointer to Dire merging objects.
Definition: Dire.h:140
Definition: DireMergingHooks.h:35
double getShowerWeight(string valKey="base")
Function to return weight of the shower evolution.
Definition: DireWeightContainer.h:156
bool initAfterBeams() override
bool Dire::init(BeamParticle* beamA, BeamParticle* beamB) {
Definition: Dire.cc:248
void registerSubObject(PhysicsBase &pb)
Register a sub object that should have its information in sync with this.
Definition: PhysicsBase.cc:56
Definition: BeamParticle.h:133
Definition: Weights.h:394
void onEndEvent(PhysicsBase::Status status) override
Definition: Dire.h:80
bool init(MergingPtr, MergingHooksPtr, PartonVertexPtr, WeightContainer *) override
Definition: Dire.h:97
Definition: DireBasics.h:374
void onBeginEvent() override
Definition: Dire.h:74
void calcWeight(double pT2, bool includeAcceptAtPT2=true, bool includeRejectAtPT2=false)
Function to calculate the weight of the shower evolution step.
Definition: DireWeightContainer.cc:381
void initShowersAndWeights()
Definition: Dire.cc:125
TimeShowerPtr getTimeShower() const override
Access the pointers to the different model components.
Definition: Dire.h:121
Definition: Dire.h:40
set< PhysicsBase * > subObjects
Definition: PhysicsBase.h:120
Status
Enumerate the different status codes the event generation can have.
Definition: PhysicsBase.h:32
void setWeightNominal(double weightNow)
The WeightContainer class.
Definition: Weights.cc:852
Container for all shower weights, including handling.
Definition: DireWeightContainer.h:82
void createPointers()
The Dire wrapper class.
Definition: Dire.cc:19
void reset()
Reset current accept/reject probabilities.
Definition: DireWeightContainer.h:115
Header for classes to set beam momentum and interaction vertex spread.
Definition: Analysis.h:20
Hooks is base class for user access to program execution.
Definition: DireHooks.h:20
Info * infoPtr
Definition: PhysicsBase.h:78
Definition: DireSplittingLibrary.h:35
Definition: ShowerModel.h:28
void setup(BeamParticle *beamA, BeamParticle *beamB)
Definition: Dire.cc:171
double weight(int i=0) const
Event weights and accumulated weight.
Definition: Info.cc:162
void initTune()
Definition: Dire.cc:54