herwig is hosted by Hepforge, IPPP Durham
Herwig  7.1.5
QTildeReconstructor.h
1 // -*- C++ -*-
2 //
3 // QTildeReconstructor.h is a part of Herwig - A multi-purpose Monte Carlo event generator
4 // Copyright (C) 2002-2017 The Herwig Collaboration
5 //
6 // Herwig is licenced under version 3 of the GPL, see COPYING for details.
7 // Please respect the MCnet academic guidelines, see GUIDELINES for details.
8 //
9 #ifndef HERWIG_QTildeReconstructor_H
10 #define HERWIG_QTildeReconstructor_H
11 //
12 // This is the declaration of the QTildeReconstructor class.
13 //
14 
15 #include "Herwig/Shower/QTilde/Base/KinematicsReconstructor.h"
16 
17 namespace Herwig {
18 
19 using namespace ThePEG;
20 
25 struct JetKinStruct {
26 
30  tShowerParticlePtr parent;
31 
35  Lorentz5Momentum p;
36 
40  Lorentz5Momentum q;
41 };
42 
46 typedef vector<JetKinStruct> JetKinVect;
47 
51 enum SystemType { UNDEFINED=-1, II, IF, F ,I };
52 
56 template<typename Value> struct ColourSinglet {
57 
58  typedef vector<ColourSinglet<Value> > VecType;
59 
60  ColourSinglet() : type(UNDEFINED) {}
61 
62  ColourSinglet(SystemType intype,Value inpart)
63  : type(intype),jets(1,inpart) {}
64 
69 
73  vector<Value> jets;
74 
75 };
76 
81 
86 
113 
114 public:
115 
119  QTildeReconstructor() : _reconopt(0), _initialBoost(0),
120  _finalStateReconOption(0),
121  _initialStateReconOption(0), _minQ(MeV) {};
122 
135  virtual bool reconstructHardJets(ShowerTreePtr hard,
136  const map<tShowerProgenitorPtr,
137  pair<Energy,double> > & pt,
138  ShowerInteraction type,
139  bool switchRecon) const;
140 
149  virtual bool reconstructDecayJets(ShowerTreePtr decay,
150  ShowerInteraction type) const;
152 
166  virtual bool deconstructDecayJets(HardTreePtr,ShowerInteraction) const;
167 
173  virtual bool deconstructHardJets(HardTreePtr,ShowerInteraction) const;
175 
176 public:
177 
184  void persistentOutput(PersistentOStream & os) const;
185 
191  void persistentInput(PersistentIStream & is, int version);
193 
200  static void Init();
201 
202 protected:
203 
221  virtual bool reconstructTimeLikeJet(const tShowerParticlePtr particleJetParent) const;
222 
231  bool reconstructSpaceLikeJet(const tShowerParticlePtr particleJetParent) const;
232 
237  bool reconstructDecayJet(const tShowerParticlePtr particleJetParent) const;
239 
249  void reconstructInitialFinalSystem(vector<ShowerProgenitorPtr>) const;
250 
255  void reconstructFinalStateSystem(bool applyBoost,
256  const LorentzRotation & toRest,
257  const LorentzRotation & fromRest,
258  vector<ShowerProgenitorPtr>) const;
259 
263  void reconstructGeneralSystem(vector<ShowerProgenitorPtr> & ShowerHardJets) const;
264 
269  void reconstructFinalFirst(vector<ShowerProgenitorPtr> & ShowerHardJets) const;
270 
275  void reconstructColourPartner(vector<ShowerProgenitorPtr> & ShowerHardJets) const;
276 
280  void reconstructColourSinglets(vector<ShowerProgenitorPtr> & ShowerHardJets,
281  ShowerInteraction type) const;
282 
287  void reconstructInitialInitialSystem(bool & applyBoost,
288  LorentzRotation & toRest,
289  LorentzRotation & fromRest,
290  vector<ShowerProgenitorPtr>) const;
292 
302  void deconstructFinalStateSystem(const LorentzRotation & toRest,
303  const LorentzRotation & fromRest,
304  HardTreePtr,
305  vector<HardBranchingPtr>,
306  ShowerInteraction) const;
307 
312  void deconstructInitialInitialSystem(bool & applyBoost,
313  LorentzRotation & toRest,
314  LorentzRotation & fromRest,
315  HardTreePtr,
316  vector<HardBranchingPtr>,
317  ShowerInteraction ) const;
318 
323  void deconstructInitialFinalSystem(HardTreePtr,
324  vector<HardBranchingPtr>,
325  ShowerInteraction ) const;
326 
327  bool deconstructGeneralSystem(HardTreePtr,
328  ShowerInteraction) const;
329 
330  bool deconstructColourSinglets(HardTreePtr,
331  ShowerInteraction) const;
332 
333  bool deconstructColourPartner(HardTreePtr,
334  ShowerInteraction) const;
336 
341  void reconstructFinalFinalOffShell(JetKinVect orderedJets, Energy2 s,
342  bool recursive) const;
343 
352  LorentzRotation solveBoost(const double k,
353  const Lorentz5Momentum & newq,
354  const Lorentz5Momentum & oldp) const;
355 
359  LorentzRotation solveBoost(const Lorentz5Momentum & newq,
360  const Lorentz5Momentum & oldq) const;
361 
365  LorentzRotation solveBoostZ(const Lorentz5Momentum & newq,
366  const Lorentz5Momentum & oldq) const;
367 
374  void boostChain(tPPtr p, const LorentzRotation & bv, tPPtr & parent) const;
375 
384  Boost solveBoostBeta( const double k, const Lorentz5Momentum & newq,
385  const Lorentz5Momentum & oldp);
386 
391  double getBeta(const double E, const double q,
392  const double Ep, const double qp) const
393  {return (q*E-qp*Ep)/(sqr(qp)+sqr(E));}
395 
411  double solveKfactor( const Energy & root_s, const JetKinVect & jets ) const;
412 
427  bool solveDecayKFactor(Energy mb,
428  const Lorentz5Momentum & n,
429  const Lorentz5Momentum & pjet,
430  const JetKinVect & jetKinematics,
431  ShowerParticlePtr partner,
432  Lorentz5Momentum ppartner[2],
433  double & k1,
434  double & k2,
435  Lorentz5Momentum & qt) const;
436 
443  double inverseRescalingFactor(vector<Lorentz5Momentum> pout,
444  vector<Energy> mon,Energy roots) const;
445 
456  bool inverseDecayRescalingFactor(vector<Lorentz5Momentum> pout,
457  vector<Energy> mon,Energy roots,
458  Lorentz5Momentum ppartner, Energy mbar,
459  double & k1, double & k2) const;
460 
467  Energy momConsEq(double k, const Energy & root_s,
468  const JetKinVect & jets) const;
469 
470 
471  void findInitialBoost(const Lorentz5Momentum & pold, const Lorentz5Momentum & pnew,
472  LorentzRotation & toRest, LorentzRotation & fromRest) const;
474 
479  template<typename Value> void findPartners(Value branch,set<Value> & done,
480  const set<Value> & branchings,
481  vector<Value> & jets) const;
482 
486  bool addIntrinsicPt(vector<ShowerProgenitorPtr>) const;
487 
496  void deepTransform(PPtr particle,const LorentzRotation & r,
497  bool match=true,PPtr original=PPtr()) const;
498 
502  Energy findMass(HardBranchingPtr) const;
503 
507  vector<double> initialStateRescaling(double x1, double x2,
508  const Lorentz5Momentum & pold,
509  const vector<Lorentz5Momentum> & p,
510  const vector<Lorentz5Momentum> & pq,
511  const vector<Energy>& highespts) const;
512 
516  vector<double> inverseInitialStateRescaling(double & x1, double & x2,
517  const Lorentz5Momentum & pold,
518  const vector<Lorentz5Momentum> & p,
519  const vector<Lorentz5Momentum> & pq) const;
520 
524  template<typename Value >
525  typename ColourSinglet<Value>::VecType identifySystems(set<Value> jets,
526  unsigned int & nnun,unsigned int & nnii,
527  unsigned int & nnif,unsigned int & nnf,
528  unsigned int & nni) const;
529 
533  template<typename Value>
534  void combineFinalState(vector<ColourSinglet<Value> > & systems) const;
535 
536 protected:
537 
544  virtual IBPtr clone() const {return new_ptr(*this);}
545 
550  virtual IBPtr fullclone() const {return new_ptr(*this);}
552 
553 protected:
554 
562  virtual void doinit();
564 
565 private:
566 
571  QTildeReconstructor & operator=(const QTildeReconstructor &) = delete;
572 
573 private:
574 
578  unsigned int _reconopt;
579 
583  unsigned int _initialBoost;
584 
589 
594 
599  Energy _minQ;
600 
604  mutable tShowerParticlePtr _progenitor;
605 
609  mutable map<tShowerProgenitorPtr,pair<Energy,double> > _intrinsic;
610 
614  mutable tShowerTreePtr _currentTree;
615 
621 
626  set<cPDPtr> _noRescale;
627 
631  mutable map<tPPtr,vector<LorentzRotation> > _boosts;
632 
636  mutable map<tShowerTreePtr,vector<LorentzRotation> > _treeBoosts;
637 };
638 
639 }
640 
641 #include "QTildeReconstructor.tcc"
642 #endif /* HERWIG_QTildeReconstructor_H */
virtual IBPtr clone() const
Make a simple clone of this object.
ThePEG::Ptr< Particle >::pointer PPtr
This class is responsible for the kinematical reconstruction after each showering step...
tShowerParticlePtr _progenitor
The progenitor of the jet currently being reconstructed.
ShowerInteraction
Handy header file to be included in all Shower classes.
map< tShowerTreePtr, vector< LorentzRotation > > _treeBoosts
Storage of the boosts applied to enable resetting after failure.
ThePEG::Ptr< InterfacedBase >::pointer IBPtr
unsigned int _finalStateReconOption
Option for the reconstruction of final state systems.
Lorentz5Momentum q
Momentum of the particle after reconstruction.
ColourSinglet< ShowerProgenitorPtr > ColourSingletSystem
Struct to store a colour singlet system of particles.
virtual IBPtr fullclone() const
Make a clone of this object, possibly modifying the cloned object to make it sane.
SystemType type
The type of system.
constexpr BinaryOpTraits< T, T >::MulT sqr(const T &x)
set< cPDPtr > _noRescale
Particles which shouldn&#39;t have their masses rescaled as set for quick access.
vector< JetKinStruct > JetKinVect
typedef for a vector of JetKinStruct
vector< Value > jets
The jets in the system.
unsigned int _initialStateReconOption
Option for the initial state reconstruction.
ThePEG::Ptr< Particle >::transient_pointer tPPtr
tShowerParticlePtr parent
Parent particle of the jet.
A simple struct to store the information we need on the showering.
tShowerTreePtr _currentTree
Current ShowerTree.
map< tPPtr, vector< LorentzRotation > > _boosts
Storage of the boosts applied to enable resetting after failure.
Lorentz5Momentum p
Momentum of the particle before reconstruction.
Energy _minQ
Minimum invariant mass for initial-final dipoles to allow the reconstruction.
map< tShowerProgenitorPtr, pair< Energy, double > > _intrinsic
Storage of the intrinsic .
double getBeta(const double E, const double q, const double Ep, const double qp) const
Compute boost parameter along z axis to get (Ep, any perp, qp) from (E, same perp, q).
-*- C++ -*-
PDVector _noRescaleVector
Particles which shouldn&#39;t have their masses rescaled as vector for the interface. ...
vector< PDPtr > PDVector
matcher< Density > match(const Generator< Density > &gen)
Indicate that the argument density should be matched to the previous one in a piecewise definition...
SystemType
Enum to identify types of colour singlet systems.
ColourSinglet< HardBranchingPtr > ColourSingletShower
Struct to store a colour singlet shower.
QTildeReconstructor()
Default constructor.
unsigned int _initialBoost
Option for the boost for initial-initial reconstruction.
unsigned int _reconopt
Option for handling the reconstruction.
This class is responsible for the kinematical reconstruction after each showering step...