herwig
is hosted by
Hepforge
,
IPPP Durham
Herwig
7.3.0
MatrixElement
FxFx
FxFxHandler.h
1
// -*- C++ -*-
2
#ifndef HERWIG_FxFxHandler_H
3
#define HERWIG_FxFxHandler_H
4
//
5
// This is the declaration of the FxFxHandler class.
6
//
7
8
#include "Herwig/Shower/QTilde/QTildeShowerHandler.h"
9
#include "Herwig/Shower/ShowerHandler.h"
10
#include "
ThePEG/Config/Pointers.h
"
11
#include "Herwig/Shower/ShowerAlpha.h"
12
#include "fastjet/PseudoJet.hh"
13
#include "fastjet/ClusterSequence.hh"
14
#include "ThePEG/Utilities/CompSelector.h"
15
#include "ThePEG/Utilities/XSecStat.h"
16
17
namespace
Herwig
{
18
class
FxFxHandler;
19
}
20
21
//declaration of thepeg ptr
22
namespace
ThePEG
{
23
ThePEG_DECLARE_POINTERS(
Herwig::FxFxHandler
,FxFxHandlerPtr);
24
}
25
26
namespace
Herwig
{
27
28
using namespace
ThePEG
;
29
30
typedef
vector< string > split_vector_type;
31
38
class
FxFxHandler
:
public
QTildeShowerHandler
{
39
43
friend
class
FxFxEventHandler
;
44
45
friend
class
FxFxReader
;
46
47
48
public
:
49
53
FxFxHandler
();
54
55
public
:
56
63
void
persistentOutput
(
PersistentOStream
& os)
const
;
64
70
void
persistentInput
(
PersistentIStream
& is,
int
version);
72
79
static
void
Init
();
80
81
protected
:
82
88
virtual
void
dofinish
();
89
95
virtual
void
doinit
();
96
101
virtual
void
doinitrun
();
103
104
public
:
109
virtual
bool
showerHardProcessVeto
()
const
;
110
114
mutable
int
npLO_
;
115
mutable
int
npNLO_;
116
120
mutable
vector<double>
ptclust_
;
121
122
protected
:
123
130
virtual
IBPtr
clone
()
const
;
131
136
virtual
IBPtr
fullclone
()
const
;
138
139
private
:
140
141
/*
142
* whether a heavy quark has been found in the merging
143
*/
144
mutable
bool
hvqfound =
false
;
145
146
/*
147
* Run MLM jet-parton matching on the 'extra' jets.
148
*/
149
bool
lightJetPartonVeto();
150
151
/*
152
* Function that calculates deltaR between a parton and a jet
153
*/
154
double
partonJetDeltaR(
ThePEG::tPPtr
partonptr, LorentzMomentum jetmom)
const
;
155
156
157
/*
158
* Function that calculates deltaR between two jets
159
*/
160
double
partonJetDeltaR(LorentzMomentum jetmom1, LorentzMomentum jetmom2)
const
;
161
165
void
getFastJets
(
double
rjet, Energy ejcut,
double
etajcut)
const
;
166
170
void
getFastJetsToMatch
(
double
rjet)
const
;
171
181
void
caldel_m
()
const
;
182
189
void
caldel_mg
()
const
;
190
200
void
caldel_hvq
()
const
;
201
205
void
getnpFxFx
()
const
;
206
210
void
getECOM
()
const
;
211
212
216
void
getptclust
()
const
;
217
222
void
erase_substr
(std::string& subject,
const
std::string& search)
const
;
223
224
229
void
getPreshowerParticles
()
const
;
230
236
void
getShoweredParticles
()
const
;
237
244
void
doSanityChecks
(
int
debugLevel)
const
;
245
250
void
getDescendents
(
PPtr
theParticle)
const
;
251
256
void
getTopRadiation
(
PPtr
theParticle)
const
;
257
262
ParticleVector
pTsort
(
ParticleVector
unsortedVec);
263
pair< vector<Energy>, vector<Lorentz5Momentum> > ETsort(vector<Energy> unsortedetjet, vector<Lorentz5Momentum> unsortedVec);
264
265
/*
266
* A function that prints a vector of Lorentz5Momenta in a fancy way
267
*/
268
void
printMomVec(vector<Lorentz5Momentum> momVec);
269
270
271
/*
272
* A probability function for varying etclus_ about the mean value
273
*/
274
Energy etclusran_(
double
petc)
const
;
275
276
private
:
277
282
static
ClassDescription<FxFxHandler>
initFxFxHandler
;
283
288
FxFxHandler
&
operator=
(
const
FxFxHandler
&) =
delete
;
289
290
private
:
291
292
297
mutable
PPair
preshowerISPs_
;
298
303
mutable
ParticleVector
preshowerFSPs_
;
304
311
mutable
ParticleVector
preshowerFSPsToDelete_
;
312
317
mutable
PPair
showeredISHs_
;
318
323
mutable
PPair
showeredISPs_
;
324
329
mutable
tPVector
showeredFSPs_
;
330
338
mutable
ParticleVector
showeredFSPsToDelete_
;
339
344
mutable
ParticleVector
partonsToMatch_
;
345
346
/*
347
* The shower progenitors
348
*/
349
350
mutable
PPtr
theProgenitor;
351
mutable
PPtr
theLastProgenitor;
352
357
mutable
tPVector
particlesToCluster_
;
358
363
mutable
PPair
showeredRems_
;
364
369
mutable
double
ECOM_
;
370
374
ShowerAlphaPtr
alphaS_
;
375
383
Energy
pdfScale_
;
384
388
Energy2
sHat_
;
389
394
double
alphaSME_
;
396
397
/*
398
* Number of rapidity segments of the calorimeter.
399
*/
400
unsigned
int
ncy_;
401
402
/*
403
* Number of phi segments of the calorimeter.
404
*/
405
unsigned
int
ncphi_;
406
407
/*
408
* Heavy flavour in WQQ,ZQQ,2Q etc (4=c, 5=b, 6=t).
409
*/
410
int
ihvy_;
411
412
/*
413
* Number of photons in the AlpGen process.
414
*/
415
int
nph_;
416
417
/*
418
* Number of higgses in the AlpGen process.
419
*/
420
int
nh_;
421
422
/*
423
* Jet ET cut to apply in jet clustering (in merging).
424
*/
425
mutable
Energy etclus_;
426
427
428
429
/*
430
* The merging mode (FxFx vs tree-level) used.
431
*/
432
int
mergemode_;
433
434
435
436
/*
437
* Allows the vetoing on heavy quark decay products to be turned off.
438
*/
439
bool
vetoHeavyQ_;
440
441
/*
442
* Allows vetoing of heavy flavour
443
*/
444
445
bool
vetoHeavyFlavour_;
446
447
448
/*
449
* Mean Jet ET cut to apply in jet clustering (in merging).
450
*/
451
Energy etclusmean_;
452
453
454
/*
455
* The jet algorithm used for parton-jet matching in the MLM procedure.
456
*/
457
int
jetAlgorithm_;
458
459
/*
460
* Allows the vetoing to be turned off completely - just for convenience.
461
*/
462
bool
vetoIsTurnedOff_;
463
464
/*
465
* Veto if there exist softer unmatched jets than matched
466
*/
467
468
bool
vetoSoftThanMatched_;
469
470
/*
471
* This flags whether the etclus_ (merging scale) should be fixed or variable according to a prob. distribution around the mean
472
*/
473
bool
etclusfixed_;
474
475
476
477
/*
478
* maximum deviation from mean Jet ET cut to apply in jet clustering (in merging).
479
*/
480
Energy epsetclus_;
481
482
483
484
/*
485
* Cone size used in jet clustering (in merging).
486
*/
487
double
rclus_;
488
489
/*
490
* Max |eta| for jets in clustering (in merging).
491
*/
492
double
etaclmax_;
493
494
/*
495
* Default 1.5 factor used to decide if a jet matches a parton
496
* in merging: if DR(parton,jet)<rclusfactor*rclus the parton
497
* and jet are said to have been matched.
498
*/
499
double
rclusfactor_;
500
501
/*
502
* Determines whether to detect the hard process or to manually determine which particles
503
* to include in the merging. If False, then the ihrd code below is used.
504
*/
505
bool
hpdetect_;
506
507
/*
508
* The AlpGen hard process code. Relation to the AlpGen process names:
509
* 1: wqq, 2: zqq, 3: wjet, 4: zjet, 5: vbjet, 6: 2Q, 8: QQh, 9: Njet,
510
* 10: wcjet, 11: phjet, 12: hjet, 13: top, 14: wphjet, 15: wphqq,
511
* 16: 2Qph.
512
*/
513
int
ihrd_;
514
515
/*
516
* The number of light jets in the AlpGen process (i.e. the 'extra' ones).
517
*/
518
int
njets_;
519
520
/*
521
* Mimimum parton-parton R-sep used for generation (used for hvq merging).
522
*/
523
double
drjmin_;
524
525
/*
526
* This flags that the highest multiplicity ME-level process is
527
* being processed.
528
*/
529
mutable
bool
highestMultiplicity_;
530
531
/*
532
* The forwards rapidity span of the calorimeter.
533
*/
534
double
ycmax_;
535
536
537
/*
538
* The backwards rapidity span of the calorimeter.
539
*/
540
double
ycmin_;
541
542
543
544
545
546
547
548
549
550
551
552
553
/*
554
* Cosine of phi values of calorimeter cell centres.
555
* Goes phi~=0 to phi~=2*pi (index = 0 ---> ncphi).
556
* ==> Cosine goes from +1 ---> +1 (index = 0 ---> ncphi).
557
*/
558
vector<double> cphcal_;
559
560
/*
561
* Sine of phi values of calorimeter cell centres.
562
* Goes phi~=0 to phi~=2*pi (index = 0 ---> ncphi).
563
* ==> Sine goes 0 -> 1 -> 0 -> -1 -> 0 (index = 0 ---> ncphi).
564
*/
565
vector<double> sphcal_;
566
567
/*
568
* Cosine of theta values of calorimeter cell centres in Y.
569
* Goes bwds th~=pi to fwds th~=0 (index = 0 ---> ncy).
570
* ==> Cosine goes from -1 ---> +1 (index = 0 ---> ncy).
571
*/
572
vector<double> cthcal_;
573
574
/*
575
* Sine of theta values of calorimeter cell centres in Y.
576
* Goes bwds th~=pi to fwds th~=0 (index = 0 ---> ncy).
577
* ==> Sine goes from 0 ---> +1 ---> 0 (index = 0 ---> ncy).
578
*/
579
vector<double> sthcal_;
580
581
/*
582
* Transverse energy deposit in a given calorimeter cell.
583
* First array index corresponds to rapidity index of cell,
584
* second array index corresponds to phi cell index.
585
*/
586
vector<vector<Energy> > et_;
587
588
/*
589
* For a given calorimeter cell this holds the index of the jet
590
* that the cell was clustered into.
591
*/
592
vector<vector<int> > jetIdx_;
593
594
/*
595
* Vector holding the Lorentz 5 momenta of each jet.
596
*/
597
mutable
vector<Lorentz5Momentum> pjet_;
598
599
/*
600
* Vector holding the Lorentz 5 momenta of each jet from ME partons
601
*/
602
mutable
vector<Lorentz5Momentum> pjetME_;
603
604
605
/*
606
* Vector holding the list of FS particles resulting from
607
* the particle input to getDescendents.
608
*/
609
mutable
ParticleVector
tmpList_;
610
611
/*
612
* Variables for the C++ translation of the calini_m(), calsim_m(),
613
* getjet_m(...) and caldel_m() functions
614
*/
615
mutable
vector<Energy> etjet_;
616
vector<Energy> etjetME_;
617
mutable
double
dely_, delphi_;
618
619
};
620
621
}
622
623
#include "ThePEG/Utilities/ClassTraits.h"
624
625
namespace
ThePEG
{
626
631
template
<>
632
struct
BaseClassTrait
<
Herwig
::FxFxHandler,1> {
634
typedef
Herwig::QTildeShowerHandler
NthBase
;
635
};
636
639
template
<>
640
struct
ClassTraits
<
Herwig
::FxFxHandler>
641
:
public
ClassTraitsBase
<Herwig::FxFxHandler> {
643
static
string
className
() {
return
"Herwig::FxFxHandler"
; }
651
static
string
library
() {
return
"HwFxFxHandler.so"
; }
652
};
653
656
}
657
658
#endif
/* HERWIG_FxFxHandler_H */
Pointers.h
Herwig::FxFxHandler
Here is the documentation of the FxFxHandler class.
Definition:
FxFxHandler.h:38
Herwig::FxFxHandler::erase_substr
void erase_substr(std::string &subject, const std::string &search) const
Erases all occurences of a substring from a string.
Herwig::FxFxHandler::showeredFSPs_
tPVector showeredFSPs_
Final-state outgoing partICLEs after shower of hard process (eventHandler()->currentEvent()->getFinal...
Definition:
FxFxHandler.h:329
Herwig::FxFxHandler::getFastJets
void getFastJets(double rjet, Energy ejcut, double etajcut) const
Find jets using the FastJet package on particlesToCluster_.
Herwig::FxFxHandler::alphaS_
ShowerAlphaPtr alphaS_
Pointer to the object calculating the strong coupling.
Definition:
FxFxHandler.h:374
Herwig::FxFxHandler::pTsort
ParticleVector pTsort(ParticleVector unsortedVec)
Sorts a given vector of particles by descending pT or ETJET.
Herwig::FxFxHandler::getptclust
void getptclust() const
get the MG5_aMC information required for tree-level merging
Herwig::FxFxHandler::doinit
virtual void doinit()
Initialize this object after the setup phase before saving an EventGenerator to disk.
Herwig::FxFxHandler::getFastJetsToMatch
void getFastJetsToMatch(double rjet) const
Find jets using the FastJet package on particlesToCluster_.
Herwig::FxFxHandler::getShoweredParticles
void getShoweredParticles() const
Get the particles from eventHandler()->currentEvent()->... filling the particle pairs showeredISHs_,...
Herwig::FxFxHandler::alphaSME_
double alphaSME_
Constant alphaS used to generate LH events - if not already using CKKW scale (ickkw = 1 in AlpGen for...
Definition:
FxFxHandler.h:394
Herwig::FxFxHandler::partonsToMatch_
ParticleVector partonsToMatch_
ONLY the final-state partons from preshowerFSPs_ that are supposed to enter the jet-parton matching.
Definition:
FxFxHandler.h:344
Herwig::FxFxHandler::getnpFxFx
void getnpFxFx() const
get the MG5_aMC information required for FxFx merging
Herwig::FxFxHandler::persistentOutput
void persistentOutput(PersistentOStream &os) const
Function used to write out object persistently.
Herwig::FxFxHandler::npLO_
int npLO_
information for FxFx merging
Definition:
FxFxHandler.h:114
Herwig::FxFxHandler::showeredISHs_
PPair showeredISHs_
Initial-state incoming hadrons after shower of hard process (eventHandler()->currentEvent()->incoming...
Definition:
FxFxHandler.h:317
Herwig::FxFxHandler::showerHardProcessVeto
virtual bool showerHardProcessVeto() const
Hook to allow vetoing of event after showering hard sub-process as in e.g.
Herwig::FxFxHandler::clone
virtual IBPtr clone() const
Make a simple clone of this object.
Herwig::FxFxHandler::particlesToCluster_
tPVector particlesToCluster_
ONLY the final-state particles from showeredFSPs_ (and maybe also showeredRems_) that are supposed to...
Definition:
FxFxHandler.h:357
Herwig::FxFxHandler::showeredISPs_
PPair showeredISPs_
Initial-state incoming partons after shower of hard process (look for partonic children of showeredIS...
Definition:
FxFxHandler.h:323
Herwig::FxFxHandler::operator=
FxFxHandler & operator=(const FxFxHandler &)=delete
The assignment operator is private and must never be called.
Herwig::FxFxHandler::FxFxHandler
FxFxHandler()
The default constructor.
Herwig::FxFxHandler::showeredRems_
PPair showeredRems_
Final-state remnants after shower of hard process (look for remnants initially in showeredFSPs_).
Definition:
FxFxHandler.h:363
Herwig::FxFxHandler::fullclone
virtual IBPtr fullclone() const
Make a clone of this object, possibly modifying the cloned object to make it sane.
Herwig::FxFxHandler::Init
static void Init()
The standard Init function used to initialize the interfaces.
Herwig::FxFxHandler::doSanityChecks
void doSanityChecks(int debugLevel) const
Allows printing of debug output and sanity checks like total momentum consrvation to be carried out.
Herwig::FxFxHandler::persistentInput
void persistentInput(PersistentIStream &is, int version)
Function used to read in object persistently.
Herwig::FxFxHandler::pdfScale_
Energy pdfScale_
Information extracted from the XComb object.
Definition:
FxFxHandler.h:383
Herwig::FxFxHandler::getECOM
void getECOM() const
get the MG5_aMC information required for FxFx merging
Herwig::FxFxHandler::initFxFxHandler
static ClassDescription< FxFxHandler > initFxFxHandler
The static object used to initialize the description of this class.
Definition:
FxFxHandler.h:282
Herwig::FxFxHandler::preshowerISPs_
PPair preshowerISPs_
Initial-state incoming partons prior to showering (i.e.
Definition:
FxFxHandler.h:297
Herwig::FxFxHandler::caldel_mg
void caldel_mg() const
Deletes particles from partonsToMatch_ and particlesToCluster_ vectors so that these contain only the...
Herwig::FxFxHandler::getDescendents
void getDescendents(PPtr theParticle) const
Given a pointer to a particle this finds all its final state descendents.
Herwig::FxFxHandler::getPreshowerParticles
void getPreshowerParticles() const
Get the particles from lastXCombPtr filling the pair preshowerISPs_ and particle pointer vector presh...
Herwig::FxFxHandler::preshowerFSPsToDelete_
ParticleVector preshowerFSPsToDelete_
Final-state outgoing partICLEs prior to showering to_be_removed from preShowerFSPs_ prior to the ligh...
Definition:
FxFxHandler.h:311
Herwig::FxFxHandler::doinitrun
virtual void doinitrun()
Initialize this object.
Herwig::FxFxHandler::caldel_hvq
void caldel_hvq() const
c++ translation of subroutine of same name from alpsho.f.
Herwig::FxFxHandler::dofinish
virtual void dofinish()
Finalize the object.
Herwig::FxFxHandler::preshowerFSPs_
ParticleVector preshowerFSPs_
Final-state outgoing partICLEs prior to showering (i.e.
Definition:
FxFxHandler.h:303
Herwig::FxFxHandler::ECOM_
double ECOM_
the COM of the incoming hadrons
Definition:
FxFxHandler.h:369
Herwig::FxFxHandler::showeredFSPsToDelete_
ParticleVector showeredFSPsToDelete_
Final-state outgoing partICLEs after shower of hard process to_be_removed from showeredFSPs_ prior to...
Definition:
FxFxHandler.h:338
Herwig::FxFxHandler::caldel_m
void caldel_m() const
Deletes particles from partonsToMatch_ and particlesToCluster_ vectors so that these contain only the...
Herwig::FxFxHandler::ptclust_
vector< double > ptclust_
information for tree-level merging
Definition:
FxFxHandler.h:120
Herwig::FxFxHandler::getTopRadiation
void getTopRadiation(PPtr theParticle) const
Accumulates all descendents of tops down to the b and W but not including them.
Herwig::FxFxHandler::sHat_
Energy2 sHat_
Centre of mass energy.
Definition:
FxFxHandler.h:388
Herwig::QTildeShowerHandler
The QTildeShowerHandler class.
Definition:
QTildeShowerHandler.h:36
ThePEG::ClassDescription
ThePEG::FxFxEventHandler
The FxFxEventHandler inherits from the general EventHandler class and administers the reading of even...
Definition:
FxFxEventHandler.h:40
ThePEG::FxFxReader
FxFxReader is an abstract base class to be used for objects which reads event files or streams from m...
Definition:
FxFxReader.h:77
ThePEG::PersistentIStream
ThePEG::PersistentOStream
Herwig
-*- C++ -*-
Definition:
BasicConsistency.h:17
ThePEG
ThePEG::tPPtr
ThePEG::Ptr< Particle >::transient_pointer tPPtr
ThePEG::PPair
pair< PPtr, PPtr > PPair
ThePEG::IBPtr
ThePEG::Ptr< InterfacedBase >::pointer IBPtr
ThePEG::tPVector
vector< tPPtr > tPVector
ThePEG::PPtr
ThePEG::Ptr< Particle >::pointer PPtr
ThePEG::ParticleVector
vector< PPtr > ParticleVector
ThePEG::BaseClassTrait
ThePEG::BaseClassTrait::NthBase
int NthBase
ThePEG::ClassTraitsBase
ThePEG::ClassTraitsBase::className
static string className()
ThePEG::ClassTraitsBase::library
static string library()
ThePEG::ClassTraits
Generated on Thu Jun 20 2024 17:50:52 for Herwig by
1.9.6