StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Dire.cc
1 // Dire.cc is a part of the PYTHIA event generator.
2 // Copyright (C) 2020 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 // Function definitions (not found in the header) for the top level
7 // Dire class.
8 
9 #include "Pythia8/Dire.h"
10 
11 namespace Pythia8 {
12 
13 //==========================================================================
14 
15 // The Dire wrapper class.
16 
17 //--------------------------------------------------------------------------
18 
19 void Dire::createPointers() {
20 
21  // Construct showers.
22  if (!weightsPtr) {
23  hasOwnWeights = true;
24  weightsPtr = new DireWeightContainer(settingsPtr);
25  }
26  if (!timesPtr) {
27  hasOwnTimes = true;
28  timesPtr = make_shared<DireTimes>( mergingHooksPtr, partonVertexPtr);
29  }
30  if (!spacePtr) {
31  hasOwnSpace = true;
32  spacePtr = make_shared<DireSpace>( mergingHooksPtr, partonVertexPtr);
33  }
34  if (!timesDecPtr) {
35  hasOwnTimesDec = true;
36  timesDecPtr = make_shared<DireTimes>( mergingHooksPtr, partonVertexPtr);
37  }
38  if (!mergingPtr) {
39  mergingPtr = make_shared<DireMerging>();
40  }
41  if (!hardProcessPtr) {
42  hasOwnHardProcess = true;
43  hardProcessPtr = new DireHardProcess();
44  }
45  if (!mergingHooksPtr) {
46  hasOwnMergingHooks = true;
47  mergingHooksPtr = make_shared<DireMergingHooks>();
48  }
49 
50 }
51 
52 //--------------------------------------------------------------------------
53 
54 void Dire::initTune() {
55 
56  initNewSettings = true;
57 
58  // Get tune id.
59  int iTune = settingsPtr->mode("Dire:Tune");
60 
61  // Default tune.
62  if (iTune == 1) {
63  // Preliminary Professor tune, dated 2017-10-10. To be used with:
64  // PDF:pSet = LHAPDF6:MMHT2014nlo68cl
65  // PDF:pHardSet = LHAPDF6:MMHT2014nlo68cl
66  // TimeShower:alphaSvalue = 0.1201
67  // SpaceShower:alphaSvalue = 0.1201
68  // ShowerPDF:usePDFalphas = on
69  // ShowerPDF:useSummedPDF = on
70  // DireSpace:forceMassiveMap = on
71  // ShowerPDF:usePDFmasses = off
72 
73  settingsPtr->readString("TimeShower:alphaSvalue = 0.1201");
74  settingsPtr->readString("SpaceShower:alphaSvalue = 0.1201");
75  settingsPtr->readString("TimeShower:alphaSorder = 2");
76  settingsPtr->readString("SpaceShower:alphaSorder = 2");
77 
78  // Tuned hadronization from e+e- data
79  settingsPtr->readString("StringPT:sigma = 0.2952");
80  settingsPtr->readString("StringZ:aLund = 0.9704");
81  settingsPtr->readString("StringZ:bLund = 1.0809");
82  settingsPtr->readString("StringZ:aExtraDiquark = 1.3490");
83  settingsPtr->readString("StringFlav:probStoUD = 0.2046");
84  settingsPtr->readString("StringZ:rFactB = 0.8321");
85  settingsPtr->readString("StringZ:aExtraSQuark = 0.0");
86  settingsPtr->readString("TimeShower:pTmin = 0.9");
87 
88  // Tuned MPI and primordial kT to LHC data (UE in dijets + Drell-Yan pT).
89  settingsPtr->readString("SpaceShower:pTmin = 0.9");
90  settingsPtr->readString("MultipartonInteractions:alphaSvalue = 0.1309");
91  settingsPtr->readString("MultipartonInteractions:pT0Ref = 1.729");
92  settingsPtr->readString("MultipartonInteractions:expPow = 1.769");
93  settingsPtr->readString("ColourReconnection:range = 2.1720");
94  settingsPtr->readString("BeamRemnants:primordialKThard = 2.2873");
95  settingsPtr->readString("BeamRemnants:primordialKTsoft = 0.25");
96  settingsPtr->readString("BeamRemnants:reducedKTatHighY = 0.47");
97 
98  }
99 
100  // For new U(1) splittings, teach Pythia new particles, if not already read
101  // from input file.
102  if ( settingsPtr->flag("TimeShower:U1newShowerByL")
103  || settingsPtr->flag("TimeShower:U1newShowerByQ")
104  || settingsPtr->flag("SpaceShower:U1newShowerByL")
105  || settingsPtr->flag("SpaceShower:U1newShowerByQ")) {
106  if (!particleDataPtr->isParticle(900032)) {
107  settingsPtr->readString("900032:all = Zp void 1 0 0 1. 0.01 0. 0. 0.");
108  settingsPtr->readString("900032:addChannel = 1 0.33 101 11 -11");
109  settingsPtr->readString("900032:addChannel = 1 0.33 101 13 -13");
110  settingsPtr->readString("900032:addChannel = 1 0.34 101 211 -211");
111  settingsPtr->readString("900032:isResonance = true");
112  }
113  if (!particleDataPtr->isParticle(900012)) {
114  settingsPtr->readString("900012:all = nup nup_bar"
115  " 1 0 0 0.0 0.0 0. 0. 0.");
116  }
117  }
118 
119  return;
120 
121 }
122 
123 //--------------------------------------------------------------------------
124 
125 void Dire::initShowersAndWeights() {
126 
127  if (isInitShower) return;
128 
129  // Construct showers.
130  if (!weightsPtr) {
131  hasOwnWeights = true;
132  weightsPtr = new DireWeightContainer(settingsPtr);
133  }
134  if (!timesPtr) {
135  hasOwnTimes = true;
136  timesPtr = make_shared<DireTimes>( mergingHooksPtr, partonVertexPtr);
137  }
138  if (!spacePtr) {
139  hasOwnSpace = true;
140  spacePtr = make_shared<DireSpace>( mergingHooksPtr, partonVertexPtr);
141  }
142  if (!timesDecPtr) {
143  hasOwnTimesDec = true;
144  timesDecPtr = make_shared<DireTimes>( mergingHooksPtr, partonVertexPtr);
145  }
146  if (!mergingPtr) {
147  mergingPtr = make_shared<DireMerging>();
148  }
149  if (!hardProcessPtr) {
150  hasOwnHardProcess = true;
151  hardProcessPtr = new DireHardProcess();
152  }
153  if (!mergingHooksPtr) {
154  hasOwnMergingHooks = true;
155  mergingHooksPtr = make_shared<DireMergingHooks>();
156  }
157 
158  mergingHooksPtr->setHardProcessPtr(hardProcessPtr);
159  mergingHooksPtr->init();
160 
161  timesPtr->setWeightContainerPtr(weightsPtr);
162  spacePtr->setWeightContainerPtr(weightsPtr);
163  timesDecPtr->setWeightContainerPtr(weightsPtr);
164 
165  isInitShower = true;
166 
167 }
168 
169 //--------------------------------------------------------------------------
170 
171 void Dire::setup(BeamParticle* beamA, BeamParticle* beamB) {
172 
173  if (isInit) return;
174 
175  // Initialise library of splitting functions.
176  if (!splittings) {
177  hasOwnSplittings = true;
178  splittings = new DireSplittingLibrary();
179  }
180 
181  // If Pythia has, for ominous reasons, not initialized the spacelike shower,
182  // retry to initialize from timelike shower beams.
183  if ( !spacePtr->isInit() && timesPtr->isInit()
184  && beamA != 0 && beamB != 0)
185  spacePtr->init( beamA, beamB );
186 
187  // Reinitialise showers to ensure that pointers are correctly set.
188  timesPtr->reinitPtr(infoPtr, mergingHooksPtr, splittings, &direInfo);
189  spacePtr->reinitPtr(infoPtr, mergingHooksPtr, splittings, &direInfo);
190  timesDecPtr->reinitPtr(infoPtr, mergingHooksPtr, splittings, &direInfo);
191 
192  // Reset Pythia masses if necessary.
193  if ( settingsPtr->flag("ShowerPDF:usePDFmasses")
194  && ( beamA != NULL || beamB != NULL) ) {
195  for (int i=1; i <= 5; ++i) {
196  // Try to get masses from the hadron beams.
197  double mPDF = (abs(beamA->id()) > 30)
198  ? beamA->mQuarkPDF(i)
199  : (abs(beamB->id()) > 30)
200  ? beamB->mQuarkPDF(i) : -1.0;
201  // If there are no hadron beams, get the masses from either beam.
202  if (beamA != NULL && mPDF < 0.)
203  mPDF = beamA->mQuarkPDF(i);
204  if (beamB != NULL && mPDF < 0.)
205  mPDF = beamB->mQuarkPDF(i);
206  if (mPDF > -1.) {
207  stringstream resetMass;
208  resetMass << i << ":m0 = " << mPDF;
209  settingsPtr->readString(resetMass.str());
210  }
211  }
212  }
213 
214  // Switch off all showering and MPI when estimating the cross section,
215  if (hooksPtr)
216  hooksPtr->initPtr( infoPtr, beamA, beamB);
217 
218  splittings->setKernelHooks(hooksPtr);
219 
220  // Initialise splitting function library here so that beam pointers
221  // are already correctly initialised.
222  splittings->init( infoPtr, beamA, beamB, &direInfo, hooksPtr);
223 
224  // Feed the splitting functions to the showers.
225  splittings->setTimesPtr(timesPtr);
226  splittings->setTimesDecPtr(timesDecPtr);
227  splittings->setSpacePtr(spacePtr);
228 
229  // Initialize splittings in showers again (!), now that splittings are
230  // properly set up.
231  timesDecPtr->initSplits();
232  timesPtr->initSplits();
233  spacePtr->initSplits();
234 
235  weightsPtr->initPtrs(beamA, beamB, settingsPtr, infoPtr, &direInfo);
236  timesDecPtr->initVariations();
237  timesPtr->initVariations();
238  spacePtr->initVariations();
239  if (mergingPtr) mergingPtr->initPtrs( weightsPtr, timesPtr,
240  spacePtr, &direInfo);
241 
242 
243 }
244 
245 //--------------------------------------------------------------------------
246 
247 //bool Dire::init(BeamParticle* beamA, BeamParticle* beamB) {
248 bool Dire::initAfterBeams() {
249 
250  if (isInit) return true;
251 
252  // Construct showers.
253  initShowersAndWeights();
254 
255  // Initialize Dire tune settings.
256  initTune();
257 
258  if ( settingsPtr->flag("Dire:doMerging")
259  || settingsPtr->flag("Dire:doMECs")
260  || settingsPtr->flag("Dire:doMEM")) {
261  settingsPtr->flag("Merging:doMerging",true);
262  settingsPtr->flag("Merging:useShowerPlugin",true);
263  }
264 
265  if ( settingsPtr->flag("Dire:doMECs")
266  || settingsPtr->flag("Dire:doMEM")) {
267  settingsPtr->parm("Merging:TMS",0.0);
268  }
269 
270  // Setup weight container (after user-defined enhance factors have been read)
271  weightsPtr->initPtrs(beamAPtr, beamBPtr, settingsPtr, infoPtr, &direInfo);
272  weightsPtr->setup();
273  setup(beamAPtr, beamBPtr);
274  isInit = true;
275 
276  printBannerSave = printBannerSave && !settingsPtr->flag("Print:quiet");
277  if (printBannerSave) printBanner();
278  printBannerSave = false;
279 
280  return isInit;
281 
282 }
283 
284 void Dire::printBanner() {
285 
286  cout << "\n"
287  << " *--------------- Welcome to the DIRE parton shower "
288  << " -------------*\n"
289  << " | "
290  << " |\n"
291  << " | Please consider citing Eur.Phys.J. C75 (2015)"
292  << " 9, 461 |\n"
293  << " | if you use this program for scientific purposes."
294  << " |\n";
295  cout << " | "
296  << " |\n";
297  cout << " *----------------------------------------"
298  << "--------------------------*" << endl;
299 
300 }
301 
302 //==========================================================================
303 
304 } // end namespace Pythia8