Line data Source code
1 : // MergingHooks.h is a part of the PYTHIA event generator.
2 : // Copyright (C) 2015 Torbjorn Sjostrand.
3 : // PYTHIA is licenced under the GNU GPL version 2, see COPYING for details.
4 : // Please respect the MCnet Guidelines, see GUIDELINES for details.
5 :
6 : // This file is written by Stefan Prestel.
7 : // Header file to allow user access to program at different stages.
8 : // HardProcess: Container class for the hard process to be merged. Holds the
9 : // bookkeeping of particles not be be reclustered
10 : // MergingHooks: Steering class for matrix element merging. Some functions can
11 : // be redefined in a derived class to have access to the merging
12 :
13 : #ifndef Pythia8_MergingHooks_H
14 : #define Pythia8_MergingHooks_H
15 :
16 : #include "Pythia8/Basics.h"
17 : #include "Pythia8/BeamParticle.h"
18 : #include "Pythia8/Event.h"
19 : #include "Pythia8/Info.h"
20 : #include "Pythia8/ParticleData.h"
21 : #include "Pythia8/PartonSystems.h"
22 : #include "Pythia8/PythiaStdlib.h"
23 : #include "Pythia8/Settings.h"
24 :
25 : namespace Pythia8 {
26 :
27 : //==========================================================================
28 :
29 : // Declaration of hard process class
30 : // This class holds information on the desired hard 2->2 process
31 : // for the merging.
32 : // This class is a container class for History class use.
33 :
34 : class HardProcess {
35 :
36 : public:
37 :
38 : // Flavour of the first incoming particle
39 : int hardIncoming1;
40 : // Flavour of the second incoming particle
41 : int hardIncoming2;
42 : // Flavours of the outgoing particles
43 : vector<int> hardOutgoing1;
44 : vector<int> hardOutgoing2;
45 : // Flavour of intermediate bosons in the hard 2->2
46 : vector<int> hardIntermediate;
47 :
48 : // Current reference event
49 : Event state;
50 : // Potential positions of outgoing particles in reference event
51 : vector<int> PosOutgoing1;
52 : vector<int> PosOutgoing2;
53 : // Potential positions of intermediate bosons in reference event
54 : vector<int> PosIntermediate;
55 :
56 : // Information on merging scale read from LHE file
57 : double tms;
58 :
59 : // Default constructor
60 0 : HardProcess(){}
61 : // Default destructor
62 0 : ~HardProcess(){}
63 :
64 : // Copy constructor
65 : HardProcess( const HardProcess& hardProcessIn )
66 : : state(hardProcessIn.state),
67 : tms(hardProcessIn.tms) {
68 : hardIncoming1 = hardProcessIn.hardIncoming1;
69 : hardIncoming2 = hardProcessIn.hardIncoming2;
70 : for(int i =0; i < int(hardProcessIn.hardOutgoing1.size());++i)
71 : hardOutgoing1.push_back( hardProcessIn.hardOutgoing1[i]);
72 : for(int i =0; i < int(hardProcessIn.hardOutgoing2.size());++i)
73 : hardOutgoing2.push_back( hardProcessIn.hardOutgoing2[i]);
74 : for(int i =0; i < int(hardProcessIn.hardIntermediate.size());++i)
75 : hardIntermediate.push_back( hardProcessIn.hardIntermediate[i]);
76 : for(int i =0; i < int(hardProcessIn.PosOutgoing1.size());++i)
77 : PosOutgoing1.push_back( hardProcessIn.PosOutgoing1[i]);
78 : for(int i =0; i < int(hardProcessIn.PosOutgoing2.size());++i)
79 : PosOutgoing2.push_back( hardProcessIn.PosOutgoing2[i]);
80 : for(int i =0; i < int(hardProcessIn.PosIntermediate.size());++i)
81 : PosIntermediate.push_back( hardProcessIn.PosIntermediate[i]);
82 : }
83 :
84 : // Constructor with path to LHE file
85 : HardProcess( string LHEfile, ParticleData* particleData) {
86 : state = Event();
87 : state.init("(hard process)", particleData);
88 : translateLHEFString(LHEfile);
89 : }
90 :
91 : // Constructor with core process input
92 : void initOnProcess( string process, ParticleData* particleData);
93 :
94 : // Constructor with path to LHE file input
95 : void initOnLHEF( string LHEfile, ParticleData* particleData);
96 :
97 : // Function to access the LHE file and read relevant information
98 : void translateLHEFString( string LHEpath);
99 :
100 : // Function to translate the process string (in MG/ME notation)
101 : void translateProcessString( string process);
102 :
103 : // Function to clear hard process information
104 : void clear();
105 :
106 : // Function to check whether the sets of candidates Pos1, Pos2, together
107 : // with the proposed candidate iPos give an allowed hard process state
108 : bool allowCandidates(int iPos, vector<int> Pos1, vector<int> Pos2,
109 : const Event& event);
110 : // Function to identify the hard subprocess in the current event
111 : void storeCandidates( const Event& event, string process);
112 : // Function to check if the particle event[iPos] matches any of
113 : // the stored outgoing particles of the hard subprocess
114 : bool matchesAnyOutgoing(int iPos, const Event& event);
115 : // Function to check if instead of the particle event[iCandidate], another
116 : // particle could serve as part of the hard process. Assumes that iCandidate
117 : // is already stored as part of the hard process.
118 : bool findOtherCandidates(int iPos, const Event& event, bool doReplace);
119 : // Function to exchange a stored hard process candidate with another choice.
120 : bool exchangeCandidates( vector<int> candidates1, vector<int> candidates2,
121 : map<int,int> further1, map<int,int> further2);
122 :
123 : // Function to get the number of coloured final state partons in the
124 : // hard process
125 : int nQuarksOut();
126 : // Function to get the number of uncoloured final state particles in the
127 : // hard process
128 : int nLeptonOut();
129 : // Function to get the number of electroweak final state bosons in the
130 : // hard process
131 : int nBosonsOut();
132 :
133 : // Function to get the number of coloured initial state partons in the
134 : // hard process
135 : int nQuarksIn();
136 : // Function to get the number of uncoloured initial state particles in the
137 : // hard process
138 : int nLeptonIn();
139 : // Function to report if a resonace decay was found in the 2->2 sub-process
140 : // of the current state
141 : int hasResInCurrent();
142 : // Function to report the number of resonace decays in the 2->2 sub-process
143 : // of the current state
144 : int nResInCurrent();
145 : // Function to report if a resonace decay was found in the 2->2 hard process
146 : bool hasResInProc();
147 : // Function to print the hard process (for debug)
148 : void list() const;
149 : // Function to print the hard process candidates in the
150 : // Matrix element state (for debug)
151 : void listCandidates() const;
152 :
153 : };
154 :
155 : //==========================================================================
156 :
157 : // MergingHooks is base class for user input to the merging procedure.
158 :
159 : class MergingHooks {
160 :
161 : public:
162 :
163 : // Constructor.
164 0 : MergingHooks() :
165 0 : doUserMergingSave(false),
166 0 : doMGMergingSave(false),
167 0 : doKTMergingSave(false),
168 0 : doPTLundMergingSave(false),
169 0 : doCutBasedMergingSave(false),
170 0 : doNL3TreeSave(false),
171 0 : doNL3LoopSave(false),
172 0 : doNL3SubtSave(false),
173 0 : doUNLOPSTreeSave(false),
174 0 : doUNLOPSLoopSave(false),
175 0 : doUNLOPSSubtSave(false),
176 0 : doUNLOPSSubtNLOSave(false),
177 0 : doUMEPSTreeSave(false),
178 0 : doUMEPSSubtSave(false),
179 0 : doEstimateXSection(false),
180 0 : doRemoveDecayProducts(false),
181 0 : doOrderHistoriesSave(true),
182 0 : doCutOnRecStateSave(false),
183 0 : doWClusteringSave(false),
184 0 : doSQCDClusteringSave(false),
185 0 : doIgnoreEmissionsSave(true),
186 0 : doIgnoreStepSave(true) {
187 0 : inputEvent = Event(); resonances.resize(0); infoPtr = 0;
188 0 : particleDataPtr = 0; partonSystemsPtr = 0;}
189 :
190 : // Make History class friend to allow access to advanced switches
191 : friend class History;
192 : // Make Pythia class friend
193 : friend class Pythia;
194 : // Make PartonLevel class friend
195 : friend class PartonLevel;
196 : // Make SpaceShower class friend
197 : friend class SpaceShower;
198 : // Make TimeShower class friend
199 : friend class TimeShower;
200 : // Make Merging class friend
201 : friend class Merging;
202 :
203 : //----------------------------------------------------------------------//
204 : // Functions that allow user interference
205 : //----------------------------------------------------------------------//
206 :
207 : // Destructor.
208 0 : virtual ~MergingHooks(){}
209 : // Function encoding the functional definition of the merging scale
210 0 : virtual double tmsDefinition( const Event& event){ return event[0].e();}
211 : // Function to dampen weights calculated from histories with lowest
212 : // multiplicity reclustered events that do not pass the ME cuts
213 : virtual double dampenIfFailCuts( const Event& inEvent ) {
214 : // Dummy statement to avoid compiler warnings
215 : if(false) cout << inEvent[0].e();
216 0 : return 1.;
217 : }
218 : // Hooks to disallow states in the construction of all histories, e.g.
219 : // because jets are below the merging scale or fail the matrix element cuts
220 : // Function to allow interference in the construction of histories
221 0 : virtual bool canCutOnRecState() { return doCutOnRecStateSave; }
222 : // Function to check reclustered state while generating all possible
223 : // histories
224 : // Function implementing check of reclustered events while constructing
225 : // all possible histories
226 : virtual bool doCutOnRecState( const Event& event ) {
227 : // Dummy statement to avoid compiler warnings.
228 : if(false) cout << event[0].e();
229 : // Count number of final state partons.
230 : int nPartons = 0;
231 0 : for( int i=0; i < int(event.size()); ++i)
232 0 : if( event[i].isFinal()
233 0 : && (event[i].isGluon() || event[i].isQuark()) )
234 0 : nPartons++;
235 : // For gg -> h, allow only histories with gluons in initial state
236 0 : if( hasEffectiveG2EW() && nPartons < 2){
237 0 : if(event[3].id() != 21 && event[4].id() != 21)
238 0 : return true;
239 : }
240 0 : return false;
241 0 : }
242 : // Function to allow not counting a trial emission.
243 0 : virtual bool canVetoTrialEmission() { return false;}
244 : // Function to check if trial emission should be rejected.
245 : virtual bool doVetoTrialEmission( const Event&, const Event& ) {
246 0 : return false; }
247 :
248 : // Function to calculate the hard process matrix element.
249 : virtual double hardProcessME( const Event& inEvent ) {
250 : // Dummy statement to avoid compiler warnings.
251 0 : if ( false ) cout << inEvent[0].e(); return 1.; }
252 :
253 : //----------------------------------------------------------------------//
254 : // Simple output functions
255 : //----------------------------------------------------------------------//
256 :
257 : // Function returning the value of the merging scale.
258 : double tms() {
259 0 : if(doCutBasedMergingSave) return 0.;
260 0 : else return tmsValueSave;
261 0 : }
262 : // Function returning the value of the Delta R_{ij} cut for
263 : // cut based merging scale definition.
264 : double dRijMS() {
265 0 : return ((tmsListSave.size() == 3) ? tmsListSave[0] : 0.);
266 : }
267 : // Function returning the value of the pT_{i} cut for
268 : // cut based merging scale definition.
269 : double pTiMS() {
270 0 : return ((tmsListSave.size() == 3) ? tmsListSave[1] : 0.);
271 : }
272 : // Function returning the value of the pT_{i} cut for
273 : // cut based merging scale definition.
274 : double QijMS() {
275 0 : return ((tmsListSave.size() == 3) ? tmsListSave[2] : 0.);
276 : }
277 : // Function returning the value of the maximal number of merged jets.
278 0 : int nMaxJets() { return nJetMaxSave;}
279 : // Function returning the value of the maximal number of merged jets,
280 : // for which NLO corrections are available.
281 0 : int nMaxJetsNLO() { return nJetMaxNLOSave;}
282 : // Function to return hard process string.
283 0 : string getProcessString() { return processSave;}
284 : // Function to return the number of outgoing partons in the core process
285 0 : int nHardOutPartons(){ return hardProcess.nQuarksOut();}
286 : // Function to return the number of outgoing leptons in the core process
287 0 : int nHardOutLeptons(){ return hardProcess.nLeptonOut();}
288 : // Function to return the number of outgoing electroweak bosons in the core
289 : // process.
290 0 : int nHardOutBosons(){ return hardProcess.nBosonsOut();}
291 : // Function to return the number of incoming partons (hadrons) in the core
292 : // process.
293 : int nHardInPartons(){ return hardProcess.nQuarksIn();}
294 : // Function to return the number of incoming leptons in the core process.
295 0 : int nHardInLeptons(){ return hardProcess.nLeptonIn();}
296 : // Function to report the number of resonace decays in the 2->2 sub-process
297 : // of the current state.
298 0 : int nResInCurrent(){ return hardProcess.nResInCurrent();}
299 : // Function to determine if user defined merging should be applied.
300 0 : bool doUserMerging(){ return doUserMergingSave;}
301 : // Function to determine if automated MG/ME merging should be applied.
302 0 : bool doMGMerging() { return doMGMergingSave;}
303 : // Function to determine if KT merging should be applied.
304 0 : bool doKTMerging() { return doKTMergingSave;}
305 : // Function to determine if PTLund merging should be applied.
306 0 : bool doPTLundMerging() { return doPTLundMergingSave;}
307 : // Function to determine if cut based merging should be applied.
308 0 : bool doCutBasedMerging() { return doCutBasedMergingSave;}
309 0 : bool doCKKWLMerging() { return (doUserMergingSave || doMGMergingSave
310 0 : || doKTMergingSave || doPTLundMergingSave || doCutBasedMergingSave); }
311 : // Functions to determine if and which part of UMEPS merging
312 : // should be applied
313 0 : bool doUMEPSTree() { return doUMEPSTreeSave;}
314 0 : bool doUMEPSSubt() { return doUMEPSSubtSave;}
315 0 : bool doUMEPSMerging() { return (doUMEPSTreeSave || doUMEPSSubtSave);}
316 : // Functions to determine if and which part of NL3 merging
317 : // should be applied
318 0 : bool doNL3Tree() { return doNL3TreeSave;}
319 : bool doNL3Loop() { return doNL3LoopSave;}
320 0 : bool doNL3Subt() { return doNL3SubtSave;}
321 0 : bool doNL3Merging() { return (doNL3TreeSave || doNL3LoopSave
322 0 : || doNL3SubtSave); }
323 : // Functions to determine if and which part of UNLOPS merging
324 : // should be applied
325 0 : bool doUNLOPSTree() { return doUNLOPSTreeSave;}
326 0 : bool doUNLOPSLoop() { return doUNLOPSLoopSave;}
327 0 : bool doUNLOPSSubt() { return doUNLOPSSubtSave;}
328 0 : bool doUNLOPSSubtNLO() { return doUNLOPSSubtNLOSave;}
329 0 : bool doUNLOPSMerging() { return (doUNLOPSTreeSave || doUNLOPSLoopSave
330 0 : || doUNLOPSSubtSave || doUNLOPSSubtNLOSave); }
331 : // Return the number clustering steps that have actually been done.
332 0 : int nRecluster() { return nReclusterSave;}
333 :
334 : //----------------------------------------------------------------------//
335 : // Output functions to analyse/prepare event for merging
336 : //----------------------------------------------------------------------//
337 :
338 : // Function to check if event contains an emission not present in the hard
339 : // process.
340 : bool isFirstEmission(const Event& event);
341 :
342 : // Function to allow effective gg -> EW boson couplings.
343 : bool hasEffectiveG2EW() {
344 0 : if ( getProcessString().compare("pp>h") == 0 ) return true;
345 0 : return false; }
346 :
347 : // Return event stripped from decay products.
348 : Event bareEvent( const Event& inputEventIn, bool storeInputEvent );
349 : // Write event with decay products attached to argument.
350 : bool reattachResonanceDecays( Event& process );
351 :
352 : // Check if particle at position iPos is part of the hard sub-system
353 : bool isInHard( int iPos, const Event& event);
354 :
355 : // Function to return the number of clustering steps for the current event
356 : int getNumberOfClusteringSteps(const Event& event);
357 :
358 : //----------------------------------------------------------------------//
359 : // Functions to steer contruction of histories
360 : //----------------------------------------------------------------------//
361 :
362 : // Function to force preferred picking of ordered histories. By default,
363 : // unordered histories will only be considered if no ordered histories
364 : // were found.
365 : void orderHistories( bool doOrderHistoriesIn) {
366 0 : doOrderHistoriesSave = doOrderHistoriesIn; }
367 : // Function to force cut on reconstructed states internally, as needed
368 : // for gg -> Higgs to ensure that e.g. uu~ -> Higgs is not constructed.
369 : void allowCutOnRecState( bool doCutOnRecStateIn) {
370 0 : doCutOnRecStateSave = doCutOnRecStateIn; }
371 :
372 : // Function to allow final state clusterings of W-bosons
373 : void doWClustering( bool doWClusteringIn ) {
374 : doWClusteringSave = doWClusteringIn; }
375 :
376 : //----------------------------------------------------------------------//
377 : // Functions used as default merging scales
378 : //----------------------------------------------------------------------//
379 :
380 : // Function to check if the input particle is a light jet, i.e. should be
381 : // checked against the merging scale defintion.
382 : bool checkAgainstCut( const Particle& particle);
383 : // Function to return the value of the merging scale function in the
384 : // current event.
385 : double tmsNow( const Event& event );
386 : // Find the minimal Lund pT between coloured partons in the event
387 : double rhoms( const Event& event, bool withColour);
388 : // Function to calculate the minimal kT in the event
389 : double kTms(const Event & event);
390 : // Find the if the event passes the Delta R_{ij}, pT_{i} and Q_{ij} cuts on
391 : // the matrix element, as needed for cut-based merging scale definition
392 : double cutbasedms( const Event& event );
393 :
394 : protected:
395 :
396 : //----------------------------------------------------------------------//
397 : // The members, switches etc.
398 : //----------------------------------------------------------------------//
399 :
400 : // Helper class doing all the core process book-keeping
401 : HardProcess hardProcess;
402 :
403 : // Pointer to various information on the generation.
404 : Info* infoPtr;
405 :
406 : // Pointer to the particle data table.
407 : ParticleData* particleDataPtr;
408 :
409 : // Pointer to the particle systems.
410 : PartonSystems* partonSystemsPtr;
411 :
412 : // AlphaS objects for alphaS reweighting use
413 : AlphaStrong AlphaS_FSRSave;
414 : AlphaStrong AlphaS_ISRSave;
415 : AlphaEM AlphaEM_FSRSave;
416 :
417 : // Saved path to LHE file for more automated merging
418 : string lheInputFile;
419 :
420 : // Flags for merging procedure definition.
421 : bool doUserMergingSave, doMGMergingSave, doKTMergingSave,
422 : doPTLundMergingSave, doCutBasedMergingSave,
423 : includeMassiveSave, enforceStrongOrderingSave, orderInRapiditySave,
424 : pickByFullPSave, pickByPoPT2Save, includeRedundantSave,
425 : pickBySumPTSave, allowColourShufflingSave, resetHardQRenSave,
426 : resetHardQFacSave;
427 : int unorderedScalePrescipSave, unorderedASscalePrescipSave,
428 : unorderedPDFscalePrescipSave, incompleteScalePrescipSave,
429 : ktTypeSave, nReclusterSave, nQuarksMergeSave;
430 : double scaleSeparationFactorSave, nonJoinedNormSave,
431 : fsrInRecNormSave, herwigAcollFSRSave, herwigAcollISRSave,
432 : pT0ISRSave, pTcutSave;
433 : bool doNL3TreeSave, doNL3LoopSave, doNL3SubtSave;
434 : bool doUNLOPSTreeSave, doUNLOPSLoopSave, doUNLOPSSubtSave,
435 : doUNLOPSSubtNLOSave;
436 : bool doUMEPSTreeSave, doUMEPSSubtSave;
437 :
438 : // Flag to only do phase space cut, rejecting events below the tms cut.
439 : bool doEstimateXSection;
440 :
441 : // Save input event in case decay products need to be detached.
442 : Event inputEvent;
443 : vector< pair<int,int> > resonances;
444 : bool doRemoveDecayProducts;
445 :
446 : // Starting scale for attaching MPI.
447 : double muMISave;
448 :
449 : // Precalculated K-factors.
450 : double kFactor0jSave;
451 : double kFactor1jSave;
452 : double kFactor2jSave;
453 :
454 : // Saved members.
455 : double tmsValueSave, DparameterSave;
456 : int nJetMaxSave;
457 : int nJetMaxNLOSave;
458 : string processSave;
459 :
460 : // List of cut values to used to define a merging scale. Ordering:
461 : // 0: DeltaR_{jet_i,jet_j,min}
462 : // 1: p_{T,jet_i,min}
463 : // 2: Q_{jet_i,jet_j,min}
464 : vector<double> tmsListSave;
465 :
466 : // INTERNAL Hooks to allow construction of all histories,
467 : // including un-ordered ones
468 : bool doOrderHistoriesSave;
469 :
470 : // INTERNAL Hooks to disallow states in the construction of all histories,
471 : // e.g. because jets are below the merging scale, of to avoid the
472 : // construction of uu~ -> Higgs histories.
473 : bool doCutOnRecStateSave;
474 :
475 : // INTERNAL Hooks to allow clustering W bosons.
476 : bool doWClusteringSave, doSQCDClusteringSave;
477 :
478 : // Store / get first scale in PDF's that Pythia should have used
479 : double muFSave;
480 : double muRSave;
481 :
482 : // Store / get factorisation scale used in matrix element calculation.
483 : double muFinMESave;
484 : double muRinMESave;
485 :
486 : // Flag to indicate trial shower usage.
487 : bool doIgnoreEmissionsSave;
488 : // Flag to indicate if events should be vetoed.
489 : bool doIgnoreStepSave;
490 : // Stored weights in case veot needs to be revoked
491 : double pTsave;
492 : double weightCKKWL1Save, weightCKKWL2Save;
493 : int nMinMPISave;
494 : // Save CKKW-L weight / O(\alpha_s) weight.
495 : double weightCKKWLSave, weightFIRSTSave;
496 :
497 : //----------------------------------------------------------------------//
498 : // Generic setup functions
499 : //----------------------------------------------------------------------//
500 :
501 : // Functions for internal use inside Pythia source code
502 : // Initialize.
503 : void init( Settings settings, Info* infoPtrIn,
504 : ParticleData* particleDataPtrIn, PartonSystems* partonSystemsPtrIn,
505 : ostream& os = cout);
506 :
507 : // Function storing candidates for the hard process in the current event
508 : // Needed in order not to cluster members of the core process
509 : void storeHardProcessCandidates(const Event& event){
510 0 : hardProcess.storeCandidates(event,getProcessString());
511 0 : }
512 : // Function to set the path to the LHE file, so that more automated merging
513 : // can be used.
514 : // Remove "_1.lhe" suffix from LHE file name.
515 : // This is done so that the HarsProcess class can access both the +0 and +1
516 : // LHE files to find both the merging scale and the core process string
517 : // Store.
518 : void setLHEInputFile( string lheFile) {
519 0 : lheInputFile = lheFile.substr(0,lheFile.size()-6); }
520 :
521 : //----------------------------------------------------------------------//
522 : // Functions for output of class members.
523 : //----------------------------------------------------------------------//
524 :
525 : // Return AlphaStrong objects
526 0 : AlphaStrong* AlphaS_FSR() { return &AlphaS_FSRSave;}
527 0 : AlphaStrong* AlphaS_ISR() { return &AlphaS_ISRSave;}
528 : AlphaEM* AlphaEM_FSR() { return &AlphaEM_FSRSave;}
529 :
530 : // Functions to return advanced merging switches
531 : // Include masses in definition of evolution pT and splitting kernels
532 0 : bool includeMassive() { return includeMassiveSave;}
533 : // Prefer strongly ordered histories
534 0 : bool enforceStrongOrdering() { return enforceStrongOrderingSave;}
535 : // Prefer histories ordered in rapidity and evolution pT
536 0 : bool orderInRapidity() { return orderInRapiditySave;}
537 : // Pick history probabilistically by full (correct) splitting probabilities
538 0 : bool pickByFull() { return pickByFullPSave;}
539 : // Pick history probabilistically, with easier form of probabilities
540 0 : bool pickByPoPT2() { return pickByPoPT2Save;}
541 : // Include redundant terms (e.g. PDF ratios) in the splitting probabilities
542 0 : bool includeRedundant(){ return includeRedundantSave;}
543 : // Pick by winner-takes-it-all, with minimum sum of scalar evolution pT
544 0 : bool pickBySumPT(){ return pickBySumPTSave;}
545 :
546 : // Prescription for combined scale of unordered emissions
547 : // 0 : use larger scale
548 : // 1 : use smaller scale
549 0 : int unorderedScalePrescip() { return unorderedScalePrescipSave;}
550 : // Prescription for combined scale used in alpha_s for unordered emissions
551 : // 0 : use combined emission scale in alpha_s weight for both (!) splittings
552 : // 1 : use original reclustered scales of each emission in alpha_s weight
553 0 : int unorderedASscalePrescip() { return unorderedASscalePrescipSave;}
554 : // Prescription for combined scale used in PDF ratios for unordered
555 : // emissions
556 : // 0 : use combined emission scale in PDFs for both (!) splittings
557 : // 1 : use original reclustered scales of each emission in PDF ratiost
558 0 : int unorderedPDFscalePrescip() { return unorderedPDFscalePrescipSave;}
559 : // Prescription for starting scale of incomplete histories
560 : // 0: use factorization scale
561 : // 1: use sHat
562 : // 2: use s
563 0 : int incompleteScalePrescip() { return incompleteScalePrescipSave;}
564 :
565 : // Allow swapping one colour index while reclustering
566 0 : bool allowColourShuffling() { return allowColourShufflingSave;}
567 :
568 : // Allow use of dynamical renormalisation scale of the core 2-> 2 process.
569 0 : bool resetHardQRen() { return resetHardQRenSave; }
570 : // Allow use of dynamical factorisation scale of the core 2-> 2 process.
571 0 : bool resetHardQFac() { return resetHardQFacSave; }
572 :
573 : // Factor by which two scales should differ to be classified strongly
574 : // ordered.
575 0 : double scaleSeparationFactor() { return scaleSeparationFactorSave;}
576 : // Absolute normalization of splitting probability for non-joined
577 : // evolution.
578 0 : double nonJoinedNorm() { return nonJoinedNormSave;}
579 : // Absolute normalization of splitting probability for final state
580 : // splittings with initial state recoiler
581 0 : double fsrInRecNorm() { return fsrInRecNormSave;}
582 : // Factor multiplying scalar evolution pT for FSR splitting, when picking
583 : // history by minimum scalar pT (see Jonathan Tully's thesis)
584 0 : double herwigAcollFSR() { return herwigAcollFSRSave;}
585 : // Factor multiplying scalar evolution pT for ISR splitting, when picking
586 : // history by minimum scalar pT (see Jonathan Tully's thesis)
587 0 : double herwigAcollISR() { return herwigAcollISRSave;}
588 : // ISR regularisation scale
589 0 : double pT0ISR() { return pT0ISRSave;}
590 : // Shower cut-off scale
591 0 : double pTcut() { return pTcutSave;}
592 :
593 : // MI starting scale.
594 0 : void muMI( double mu) { muMISave = mu; }
595 0 : double muMI() { return muMISave;}
596 :
597 : // Full k-Factor for current event
598 : double kFactor(int njet = 0) {
599 0 : return (njet == 0) ? kFactor0jSave
600 0 : :(njet == 1) ? kFactor1jSave
601 0 : : kFactor2jSave;
602 : }
603 : // O(\alhpa_s)-term of the k-Factor for current event
604 : double k1Factor( int njet = 0) {
605 0 : return (kFactor(njet) - 1)/infoPtr->alphaS();
606 : }
607 :
608 : // Function to return if construction of histories is biased towards ordered
609 : // histories.
610 0 : bool orderHistories() { return doOrderHistoriesSave;}
611 :
612 : // INTERNAL Hooks to disallow states in the construction of all histories,
613 : // e.g. because jets are below the merging scale, of to avoid the
614 : // construction of uu~ -> Higgs histories.
615 0 : bool allowCutOnRecState() { return doCutOnRecStateSave;}
616 :
617 : // INTERNAL Hooks to allow clustering W bosons.
618 0 : bool doWClustering() { return doWClusteringSave;}
619 : // INTERNAL Hooks to allow clustering clustering of gluons to squarks.
620 0 : bool doSQCDClustering() { return doSQCDClusteringSave;}
621 :
622 : // Store / get first scale in PDF's that Pythia should have used
623 0 : double muF() { return (muFSave > 0.) ? muFSave : infoPtr->QFac();}
624 0 : double muR() { return (muRSave > 0.) ? muRSave : infoPtr->QRen();}
625 : // Store / get factorisation scale used in matrix element calculation.
626 : double muFinME() {
627 0 : double mu = atof((char*)infoPtr->getEventAttribute("muf2",true).c_str());
628 0 : return (muFinMESave > 0.) ? muFinMESave
629 0 : : (mu > 0.) ? sqrt(mu) : infoPtr->QFac();
630 0 : }
631 : double muRinME() {
632 0 : double mu = atof((char*)infoPtr->getEventAttribute("mur2",true).c_str());
633 0 : return (muRinMESave > 0.) ? muRinMESave
634 0 : : (mu > 0.) ? sqrt(mu) : infoPtr->QRen();
635 0 : }
636 :
637 : //----------------------------------------------------------------------//
638 : // Functions to steer shower evolution
639 : //----------------------------------------------------------------------//
640 :
641 : // Flag to indicate trial shower usage.
642 : void doIgnoreEmissions( bool doIgnoreIn ) {
643 0 : doIgnoreEmissionsSave = doIgnoreIn;
644 0 : }
645 : // Function to allow not counting a trial emission.
646 0 : bool canVetoEmission() { return !doIgnoreEmissionsSave; }
647 : // Function to check if emission should be rejected.
648 : bool doVetoEmission( const Event& );
649 :
650 : // Flag to indicate if events should be vetoed.
651 0 : void doIgnoreStep( bool doIgnoreIn ) { doIgnoreStepSave = doIgnoreIn; }
652 : // Function to allow event veto.
653 : bool canVetoStep() { return !doIgnoreStepSave; }
654 : // Function to check event veto.
655 : bool doVetoStep( const Event& process, const Event& event,
656 : bool doResonance = false );
657 :
658 : // Stored weights in case veot needs to be revoked
659 0 : void storeWeights( double weight ){ weightCKKWL1Save = weightCKKWL2Save
660 0 : = weight; }
661 :
662 : // Set starting scales
663 : bool setShowerStartingScales( bool isTrial, bool doMergeFirstEmm,
664 : double& pTscaleIn, const Event& event,
665 : double& pTmaxFSRIn, bool& limitPTmaxFSRin,
666 : double& pTmaxISRIn, bool& limitPTmaxISRin,
667 : double& pTmaxMPIIn, bool& limitPTmaxMPIin );
668 0 : void nMinMPI( int nMinMPIIn ) { nMinMPISave = nMinMPIIn; }
669 0 : int nMinMPI() { return nMinMPISave;}
670 :
671 : //----------------------------------------------------------------------//
672 : // Functions for internal merging scale definions
673 : //----------------------------------------------------------------------//
674 :
675 : // Function to calculate the kT separation between two particles
676 : double kTdurham(const Particle& RadAfterBranch,
677 : const Particle& EmtAfterBranch, int Type, double D );
678 : // Function to compute "pythia pT separation" from Particle input
679 : double rhoPythia(const Particle& RadAfterBranch,
680 : const Particle& EmtAfterBranch, const Particle& RecAfterBranch,
681 : int ShowerType);
682 : // Function to find a colour (anticolour) index in the input event,
683 : // used to find colour-connected recoilers
684 : int findColour(int col, int iExclude1, int iExclude2,
685 : const Event& event, int type, bool isHardIn);
686 : // Function to compute Delta R separation from 4-vector input
687 : double deltaRij(Vec4 jet1, Vec4 jet2);
688 :
689 : //----------------------------------------------------------------------//
690 : // Functions for weight management
691 : //----------------------------------------------------------------------//
692 :
693 : // Function to get the CKKW-L weight for the current event
694 : double getWeightNLO() { return (weightCKKWLSave - weightFIRSTSave);}
695 : // Return CKKW-L weight.
696 0 : double getWeightCKKWL() { return weightCKKWLSave; }
697 : // Return O(\alpha_s) weight.
698 : double getWeightFIRST() { return weightFIRSTSave; }
699 : // Set CKKW-L weight.
700 : void setWeightCKKWL( double weightIn){
701 0 : weightCKKWLSave = weightIn;
702 0 : infoPtr->setWeightCKKWL(weightIn); }
703 : // Set O(\alpha_s) weight.
704 : void setWeightFIRST( double weightIn){
705 0 : weightFIRSTSave = weightIn;
706 0 : infoPtr->setWeightFIRST(weightIn); }
707 :
708 : };
709 :
710 : //==========================================================================
711 :
712 : } // end namespace Pythia8
713 :
714 : #endif // Pythia8_MergingHooks_H
|
