66#include "Generators/GeneratorPythia8.h"
77#include "fairlogger/Logger.h"
88#include "TRandom3.h"
9+ #include "TParticlePDG.h"
10+ #include "TDatabasePDG.h"
911#include "TSystem.h"
10- #include <algorithm>
12+ #include "TMath.h"
1113#include <cmath>
12- #include <cstdlib>
13- #include <cstdint>
1414#include <vector>
15+ #include <fstream>
1516#include <string>
17+ using namespace Pythia8 ;
1618#endif
1719
1820R__ADD_INCLUDE_PATH ($O2DPG_MC_CONFIG_ROOT )
1921#include "MC/config/common/external/generator/CoalescencePythia8.h"
22+ /// First version of the simple coalescence generator based PYTHIA8
2023
21- class GeneratorPythia8Coalescence
22- : public o2 ::eventgen ::GeneratorPythia8
24+ class GeneratorPythia8Coalescence : public o2 ::eventgen ::GeneratorPythia8
2325{
2426public :
25- GeneratorPythia8Coalescence (
26- std ::vector < int > pdgList ,
27- int input_trigger_ratio = 1 ,
28- double coal_momentum = 0.4 ,
29- bool trivialCoal = false)
30- : o2 ::eventgen ::GeneratorPythia8 (),
31- mRequestedPdgList (std ::move (pdgList )),
32- mInverseTriggerRatio (input_trigger_ratio ),
33- mCoalMomentum (coal_momentum ),
34- mTrivialCoal (trivialCoal )
27+ /// Constructor
28+ GeneratorPythia8Coalescence (std ::vector < unsigned int > pdgList , int input_trigger_ratio = 1 , double coal_momentum = 0.4 )
29+ : o2 ::eventgen ::GeneratorPythia8 ()
3530 {
36- if (mRequestedPdgList .empty () ||
37- mInverseTriggerRatio < 1 ) {
38- LOG (fatal ) << "Invalid coalescence configuration" ;
39- }
40-
41- // The afterburner expects positive PDG codes.
42- // Preserve the signed list for event acceptance.
43- for (auto pdg : mRequestedPdgList ) {
44- unsigned int absPdg =
45- static_cast < unsigned int > (std ::abs (pdg ));
46-
47- if (std ::find (mPdgList .begin (),
48- mPdgList .end (),
49- absPdg ) == mPdgList .end ()) {
50- mPdgList .push_back (absPdg );
51- }
52- }
53-
54- LOG (info ) << "Coalescence generator: p0 = "
55- << mCoalMomentum
56- << ", trivialCoal = " << mTrivialCoal
57- << ", trigger ratio = " << mInverseTriggerRatio ;
31+ fmt ::printf (">> Coalescence generator %d\n" , input_trigger_ratio );
32+ mInverseTriggerRatio = input_trigger_ratio ;
33+ mCoalMomentum = coal_momentum ;
34+ mPdgList = pdgList ;
5835 }
59-
60- ~GeneratorPythia8Coalescence () override = default ;
36+ /// Destructor
37+ ~GeneratorPythia8Coalescence () = default ;
6138
6239 bool Init () override
6340 {
@@ -69,108 +46,55 @@ public:
6946protected :
7047 bool generateEvent () override
7148 {
72- if (mGeneratedEvents % mInverseTriggerRatio == 0 ) {
73-
49+ fmt ::printf (">> Generating event %d\n" , mGeneratedEvents );
50+ // Simple straightforward check to alternate generators
51+ if (mGeneratedEvents % mInverseTriggerRatio == 0 )
52+ {
53+ fmt ::printf (">> Generating coalescence event %d\n" , mGeneratedEvents );
7454 bool genOk = false;
75- uint64_t localCounter = 0 ;
76-
77- while (! genOk ) {
78-
79- if (! GeneratorPythia8 :: generateEvent ()) {
80- continue ;
55+ int localCounter { 0 } ;
56+ while (! genOk )
57+ {
58+ if ( GeneratorPythia8 :: generateEvent ())
59+ {
60+ genOk = CoalescencePythia8 ( mPythia . event , mPdgList , false, mCoalMomentum ) ;
8161 }
82-
83- ++ localCounter ;
84-
85- bool coalHappened = CoalescencePythia8 (
86- mPythia .event ,
87- mPdgList ,
88- mTrivialCoal ,
89- mCoalMomentum
90- );
91-
92- if (!coalHappened ) {
93- continue ;
94- }
95-
96- // Accept only events containing a requested
97- // nucleus or antinucleus.
98- genOk = containsRequestedNucleus (mPythia .event );
62+ localCounter ++ ;
9963 }
100-
101- LOG (info ) << "Coalescence successful after "
102- << localCounter << " generations" ;
103-
64+ fmt ::printf (">> Coalescence successful after %i generations\n" , localCounter );
65+ std ::cout << std ::endl << std ::endl ;
10466 notifySubGenerator (1 );
105-
106- } else {
107-
67+ }
68+ else
69+ {
70+ fmt ::printf (">> Generating minimum-bias event %d\n" , mGeneratedEvents );
71+ // Generate minimum-bias event
10872 bool genOk = false;
109-
110- while (! genOk ) {
73+ while (! genOk )
74+ {
11175 genOk = GeneratorPythia8 ::generateEvent ();
11276 }
113-
11477 notifySubGenerator (0 );
11578 }
116-
117- ++ mGeneratedEvents ;
79+ mGeneratedEvents ++ ;
11880 return true;
11981 }
12082
121- bool containsRequestedNucleus (
122- const Pythia8 ::Event & event ) const
123- {
124- for (int i = 0 ; i < event .size (); ++ i ) {
125-
126- if (!event [i ].isFinal ()) {
127- continue ;
128- }
129-
130- const int pdg = event [i ].id ();
131-
132- if (std ::find (mRequestedPdgList .begin (),
133- mRequestedPdgList .end (),
134- pdg ) != mRequestedPdgList .end ()) {
135- return true;
136- }
137- }
138-
139- return false;
140- }
14183
14284private :
143- // Signed PDG codes requested by the user
144- std ::vector < int > mRequestedPdgList ;
145-
146- // Absolute PDG codes used by the afterburner
147- std ::vector < unsigned int > mPdgList ;
148-
149- int mInverseTriggerRatio = 1 ;
150- double mCoalMomentum = 0.4 ;
151- bool mTrivialCoal = false;
152-
153- uint64_t mGeneratedEvents = 0 ;
85+ std ::vector < unsigned int > mPdgList ; /// list of pdg codes to be generated
86+ float mCoalMomentum = 0.4 ; /// coalescence momentum
87+ uint64_t mGeneratedEvents = 0 ; /// number of events generated so far
88+ int mInverseTriggerRatio = 1 ; /// injection gap
15489};
15590
156- /// Generator factory
157- FairGenerator * generateCoalescence (
158- std ::vector < int > pdgList ,
159- int input_trigger_ratio = 1 ,
160- double coal_momentum = 0.4 ,
161- bool trivialCoal = false)
91+ ///___________________________________________________________
92+ FairGenerator * generateCoalescence (std ::vector < unsigned int > pdgList , int input_trigger_ratio , double coal_momentum = 0.4 )
16293{
163- auto myGen = new GeneratorPythia8Coalescence (
164- pdgList ,
165- input_trigger_ratio ,
166- coal_momentum ,
167- trivialCoal
168- );
169-
94+ auto myGen = new GeneratorPythia8Coalescence (pdgList , input_trigger_ratio , coal_momentum );
17095 auto seed = (gRandom -> TRandom ::GetSeed () % 900000000 );
171-
17296 myGen -> readString ("Random:setSeed on" );
17397 myGen -> readString ("Random:seed " + std ::to_string (seed ));
174-
17598 return myGen ;
17699}
100+ ///___________________________________________________________
0 commit comments