KNOPythiaHadronization.cxx
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*
3  Copyright (c) 2003-2019, The GENIE Collaboration
4  For the full text of the license visit http://copyright.genie-mc.org
5  or see $GENIE/LICENSE
6 
7  Author: Costas Andreopoulos <costas.andreopoulos \at stfc.ac.uk>
8  University of Liverpool & STFC Rutherford Appleton Lab
9 
10  For the class documentation see the corresponding header file.
11 
12  Important revisions after version 2.0.0 :
13  @ Feb 10, 2011
14  Fixed a bug reported by Torben Ferber affecting the KNO -> PYTHIA
15  model transition (the order was reversed!)
16 
17 */
18 //____________________________________________________________________________
19 
20 #include <cstdlib>
21 
22 #include <TLorentzVector.h>
23 #include <TClonesArray.h>
24 #include <TH1D.h>
25 
33 
34 using namespace genie;
35 using namespace genie::constants;
36 
37 //____________________________________________________________________________
39 HadronizationModelI("genie::KNOPythiaHadronization")
40 {
41 
42 }
43 //____________________________________________________________________________
45 HadronizationModelI("genie::KNOPythiaHadronization", config)
46 {
47 
48 }
49 //____________________________________________________________________________
51 {
52 
53 }
54 //____________________________________________________________________________
56 {
57 
58 }
59 //____________________________________________________________________________
61  const Interaction * interaction) const
62 {
63 // Generate the hadronic system using either the KNO-based or PYTHIA/JETSET
64 // hadronization models according to the specified transition scheme
65 
66  double W = interaction->Kine().W();
67  LOG("HybridHad", pINFO) << "W = " << W << " GeV";
68 
69  if(W <= kNucleonMass+kPionMass) {
70  LOG("HybridHad", pWARN)
71  << "Low invariant mass, W = " << W << " GeV! Returning a null list";
72  return 0;
73  }
74 
75  //-- Init event weight (to be set if producing weighted events)
76  fWeight = 1.;
77 
78  //-- Select hadronizer
79  const HadronizationModelI * hadronizer = this->SelectHadronizer(interaction);
80 
81  //-- Run the selected hadronizer
82  TClonesArray * particle_list = hadronizer->Hadronize(interaction);
83 
84  //-- Update the weight
85  fWeight = hadronizer->Weight();
86 
87  return particle_list;
88 }
89 //____________________________________________________________________________
91  const Interaction * interaction) const
92 {
93  //-- Select hadronizer
94  const HadronizationModelI * hadronizer = this->SelectHadronizer(interaction);
95 
96  //-- Run the selected hadronizer
97  PDGCodeList * pdgv = hadronizer->SelectParticles(interaction);
98 
99  return pdgv;
100 }
101 //____________________________________________________________________________
103  const Interaction * interaction, Option_t * opt) const
104 {
105  //-- Select hadronizer
106  const HadronizationModelI * hadronizer = this->SelectHadronizer(interaction);
107 
108  //-- Run the selected hadronizer
109  TH1D * mprob = hadronizer->MultiplicityProb(interaction,opt);
110 
111  return mprob;
112 }
113 //____________________________________________________________________________
115 {
116  return fWeight;
117 }
118 //____________________________________________________________________________
120  const Interaction * interaction) const
121 {
122  const HadronizationModelI * hadronizer = 0;
124  double W = 0;
125 
126  switch(fMethod) {
127 
128  // ** KNO-only
129  case(0) :
130  hadronizer = fKNOHadronizer;
131  break;
132 
133  // ** PYTHIA/JETSET-only
134  case(1) :
135  hadronizer = fPythiaHadronizer;
136  break;
137 
138  // ** KNO-only @ W < Wmin
139  // ** PYTHIA/JETSET-only @ W > Wmax
140  // ** Smooth linear transition in [Wmin,Wmax]
141  case(2) :
142  W = interaction->Kine().W();
143  if (W <= fWminTrWindow) hadronizer = fKNOHadronizer;
144  else if (W > fWmaxTrWindow) hadronizer = fPythiaHadronizer;
145  else {
146  // Transition window
147  double R = rnd->RndHadro().Rndm();
149  if(R>f) hadronizer = fKNOHadronizer;
150  else hadronizer = fPythiaHadronizer;
151  }
152  break;
153 
154  default :
155  LOG("HybridHad", pFATAL)
156  << "Unspecified transition method: " << fMethod;
157  exit(1);
158  }
159 
160  if(!hadronizer) {
161  LOG("HybridHad", pFATAL) << "Null hadronizer!!";
162  exit(1);
163  }
164 
165  LOG("HybridHad", pINFO) << "Selected hadronizer: " << hadronizer->Id();
166  return hadronizer;
167 }
168 //____________________________________________________________________________
170 {
171  Algorithm::Configure(config);
172  this->LoadConfig();
173 }
174 //____________________________________________________________________________
176 {
177  Algorithm::Configure(config);
178  this->LoadConfig();
179 }
180 //____________________________________________________________________________
182 {
183 // Read configuration options or set defaults
184 
185  // Load the requested hadronizers
186  fKNOHadronizer =
187  dynamic_cast<const HadronizationModelI *> (
188  this->SubAlg("KNO-Hadronizer"));
190  dynamic_cast<const HadronizationModelI *> (
191  this->SubAlg("PYTHIA-Hadronizer"));
192 
194 
195  // Get transition method
196  fMethod = 2 ;
197  GetParam( "TransMethod", fMethod, false ) ;
198 
199 
200  // Get transition scheme specific config
201  if(fMethod==2) {
202 
203  GetParam( "KNO2PYTHIA-Wmin", fWminTrWindow ) ;
204 
205  GetParam( "KNO2PYTHIA-Wmax", fWmaxTrWindow ) ;
206 
207  }
208 }
209 //____________________________________________________________________________
double W(bool selected=false) const
Definition: Kinematics.cxx:167
Basic constants.
int fMethod
KNO -> PYTHIA transition method.
THE MAIN GENIE PROJECT NAMESPACE
Definition: GeneratorBase.h:8
static const double kNucleonMass
Definition: Constants.h:78
static RandomGen * Instance()
Access instance.
Definition: RandomGen.cxx:79
void Initialize(void) const
define the HadronizationModelI interface
#define pFATAL
Definition: Messenger.h:57
const HadronizationModelI * SelectHadronizer(const Interaction *) const
TClonesArray * Hadronize(const Interaction *) const
Definition: config.py:1
void Configure(const Registry &config)
PDGCodeList * SelectParticles(const Interaction *) const
A singleton holding random number generator classes. All random number generation in GENIE should tak...
Definition: RandomGen.h:30
A list of PDG codes.
Definition: PDGCodeList.h:33
double fWminTrWindow
min W in transition region (pure KNO < Wmin)
double fWeight
weight for generated event
Summary information for an interaction.
Definition: Interaction.h:56
double fWmaxTrWindow
max W in transition region (pure PYTHIA > Wmax)
virtual PDGCodeList * SelectParticles(const Interaction *) const =0
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:97
TH1D * MultiplicityProb(const Interaction *, Option_t *opt="") const
const HadronizationModelI * fPythiaHadronizer
PYTHIA Hadronizer.
const Kinematics & Kine(void) const
Definition: Interaction.h:71
virtual void Configure(const Registry &config)
Definition: Algorithm.cxx:70
#define R(x)
virtual TClonesArray * Hadronize(const Interaction *) const =0
#define pINFO
Definition: Messenger.h:63
virtual TH1D * MultiplicityProb(const Interaction *, Option_t *opt="") const =0
#define pWARN
Definition: Messenger.h:61
const HadronizationModelI * fKNOHadronizer
KNO Hadronizer.
TRandom3 & RndHadro(void) const
rnd number generator used by hadronization models
Definition: RandomGen.h:54
virtual const AlgId & Id(void) const
Get algorithm ID.
Definition: Algorithm.h:98
static const double kPionMass
Definition: Constants.h:74
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:66
Pure abstract base class. Defines the HadronizationModelI interface to be implemented by any algorith...
exit(0)
assert(nhit_max >=nhit_nbins)
bool GetParam(const RgKey &name, T &p, bool is_top_call=true) const
#define W(x)
virtual double Weight(void) const =0
const Algorithm * SubAlg(const RgKey &registry_key) const
Definition: Algorithm.cxx:353