ThePEG  1.8.0
Cuts.h
1 // -*- C++ -*-
2 //
3 // Cuts.h is a part of ThePEG - Toolkit for HEP Event Generation
4 // Copyright (C) 1999-2011 Leif Lonnblad
5 //
6 // ThePEG is licenced under version 2 of the GPL, see COPYING for details.
7 // Please respect the MCnet academic guidelines, see GUIDELINES for details.
8 //
9 #ifndef THEPEG_Cuts_H
10 #define THEPEG_Cuts_H
11 //
12 // This is the declaration of the Cuts class.
13 //
14 
15 #include "ThePEG/Interface/Interfaced.h"
16 #include "Cuts.fh"
17 #include "OneCutBase.h"
18 #include "TwoCutBase.h"
19 #include "MultiCutBase.h"
20 #include "JetFinder.h"
21 
22 namespace ThePEG {
23 
52 class Cuts: public Interfaced {
53 
54 public:
55 
59  typedef vector<OneCutPtr> OneCutVector;
60 
64  typedef vector<TwoCutPtr> TwoCutVector;
65 
69  typedef vector<MultiCutPtr> MultiCutVector;
70 
71 public:
72 
78  Cuts(Energy MhatMin=2*GeV);
79 
83  virtual ~Cuts();
85 
86 public:
87 
98  virtual void initialize(Energy2 smax, double Y);
99 
106  virtual void initEvent();
107 
123  virtual bool
124  initSubProcess(Energy2 shat, double yhat, bool mirror = false) const;
126 
136  virtual bool passCuts(const tcPDVector & ptype, const vector<LorentzMomentum> & p,
137  tcPDPtr t1 = tcPDPtr(), tcPDPtr t2 = tcPDPtr()) const;
138 
146  bool passCuts(const tcPVector & p,
147  tcPDPtr t1 = tcPDPtr(), tcPDPtr t2 = tcPDPtr()) const;
148 
154  bool passCuts(const SubProcess & sub) const;
155 
160  bool passCuts(const Collision & coll) const;
162 
172  Energy2 minSij(tcPDPtr pi, tcPDPtr pj) const;
173 
182  Energy2 minTij(tcPDPtr pi, tcPDPtr po) const;
183 
191  double minDeltaR(tcPDPtr pi, tcPDPtr pj) const;
192 
205  Energy minKTClus(tcPDPtr pi, tcPDPtr pj) const;
206 
214  double minDurham(tcPDPtr pi, tcPDPtr pj) const;
215 
222  Energy minKT(tcPDPtr p) const;
223 
230  double minEta(tcPDPtr p) const;
231 
238  double maxEta(tcPDPtr p) const;
239 
246  double minRapidityMax(tcPDPtr p) const;
247 
254  double maxRapidityMin(tcPDPtr p) const;
255 
261  double minYStar(tcPDPtr p) const;
262 
268  double maxYStar(tcPDPtr p) const;
269 
277  Energy2 minS(const tcPDVector & pv) const;
278 
286  Energy2 maxS(const tcPDVector & pv) const;
288 
295  template <typename T>
296  vector<typename Ptr<T>::transient_const_pointer>
297  oneCutObjects() const;
298 
303  template <typename T>
304  vector<typename Ptr<T>::transient_const_pointer>
305  twoCutObjects() const;
306 
311  template <typename T>
312  vector<typename Ptr<T>::transient_const_pointer>
313  multiCutObjects() const;
314 
319  const OneCutVector& oneCuts() const { return theOneCuts; }
320 
325  const TwoCutVector& twoCuts() const { return theTwoCuts; }
326 
331  const MultiCutVector& multiCuts() const { return theMultiCuts; }
332 
337 
341  void add(tOneCutPtr c) { theOneCuts.push_back(c); }
342 
346  void add(tTwoCutPtr c) { theTwoCuts.push_back(c); }
347 
351  void add(tMultiCutPtr c) { theMultiCuts.push_back(c); }
353 
354 public:
355 
362  Energy2 SMax() const { return theSMax; }
363 
364 
369  double Y() const { return theY; }
370 
375  Energy2 currentSHat() const { return theCurrentSHat; }
376 
382  double currentYHat() const { return theCurrentYHat; }
383 
385 
391  Energy2 sHatMin() const { return max(sqr(theMHatMin), theX1Min*theX2Min*SMax()); }
392 
396  Energy2 sHatMax() const { return min(sqr(theMHatMax), theX1Max*theX2Max*SMax()); }
397 
401  bool sHat(Energy2 sh) const {
402  return sh > sHatMin() && sh <= sHatMax()*(1.0 + 1000.0*Constants::epsilon);
403  }
404 
408  Energy mHatMin() const { return max(theMHatMin, sqrt(theX1Min*theX2Min*SMax())); }
409 
413  Energy mHatMax() const { return min(theMHatMax, sqrt(theX1Max*theX2Max*SMax())); }
414 
419  double yHatMin() const;
420 
425  double yHatMax() const;
426 
430  bool yHat(double y) const;
431 
436  double x1Min() const;
437 
442  double x1Max() const;
443 
447  bool x1(double x) const;
448 
453  double x2Min() const;
454 
459  double x2Max() const;
460 
464  bool x2(double x) const;
465 
470  Energy2 scaleMin() const { return theScaleMin; }
471 
476  Energy2 scaleMax() const { return theScaleMax; }
477 
481  bool scale(Energy2 Q2) const { return Q2 > scaleMin() && Q2 < scaleMax(); }
482 
487  bool subMirror() const { return theSubMirror; }
489 
490 public:
491 
495  virtual void describe() const;
496 
497 protected:
498 
505  virtual void doinitrun();
507 
508 public:
509 
516  void persistentOutput(PersistentOStream & os) const;
517 
523  void persistentInput(PersistentIStream & is, int version);
525 
532  static void Init();
533 
534 protected:
535 
542  virtual IBPtr clone() const;
543 
548  virtual IBPtr fullclone() const;
550 
551 private:
552 
556  Energy maxMHatMin() const;
557 
561  Energy minMHatMax() const;
562 
566  double maxYHatMin() const;
567 
571  double minYHatMax() const;
572 
576  double maxX1Min() const;
577 
581  double minX1Max() const;
582 
586  double maxX2Min() const;
587 
591  double minX2Max() const;
592 
596  Energy2 maxScaleMin() const;
597 
601  Energy2 minScaleMax() const;
602 
603 private:
604 
609  Energy2 theSMax;
610 
615  double theY;
616 
621  mutable Energy2 theCurrentSHat;
622 
628  mutable double theCurrentYHat;
629 
634 
639 
644  double theYHatMin;
645 
650  double theYHatMax;
651 
656  double theX1Min;
657 
662  double theX1Max;
663 
668  double theX2Min;
669 
674  double theX2Max;
675 
680  Energy2 theScaleMin;
681 
686  Energy2 theScaleMax;
687 
692  OneCutVector theOneCuts;
693 
698  TwoCutVector theTwoCuts;
699 
704  MultiCutVector theMultiCuts;
705 
711 
716  mutable bool theSubMirror;
717 
718 private:
719 
725 
730  Cuts & operator=(const Cuts &);
731 
732 };
733 
734 }
735 
736 #include "ThePEG/Utilities/ClassTraits.h"
737 
738 namespace ThePEG {
739 
744 template <>
745 struct BaseClassTrait<Cuts,1> {
747  typedef Interfaced NthBase;
748 };
749 
752 template <>
753 struct ClassTraits<Cuts>
754  : public ClassTraitsBase<Cuts> {
756  static string className() { return "ThePEG::Cuts"; }
757 };
758 
761 }
762 
763 #endif /* THEPEG_Cuts_H */
Energy maxMHatMin() const
Helper function used by the interface.
Energy2 currentSHat() const
The invariant mass squared of the hard sub-process of the event being considered. ...
Definition: Cuts.h:375
Energy2 minSij(tcPDPtr pi, tcPDPtr pj) const
Return the minimum allowed squared invariant mass of two outgoing partons of type pi and pj...
PersistentIStream is used to read persistent objects from a stream where they were previously written...
void persistentInput(PersistentIStream &is, int version)
Function used to read in object persistently.
bool yHat(double y) const
Check if the given is within the cuts.
double currentYHat() const
The total rapidity of hard sub-process (wrt.
Definition: Cuts.h:382
double minYHatMax() const
Helper function used by the interface.
double theY
The total rapidity of the colliding particles corresponding to the maximum invariant mass squared...
Definition: Cuts.h:615
double maxEta(tcPDPtr p) const
Return the maximum allowed pseudo-rapidity of an outgoing parton of the given type.
double Y() const
The total rapidity of the colliding particles corresponding to the maximum invariant mass squared...
Definition: Cuts.h:369
virtual bool passCuts(const tcPDVector &ptype, const vector< LorentzMomentum > &p, tcPDPtr t1=tcPDPtr(), tcPDPtr t2=tcPDPtr()) const
Check if the outgoing particles, with the given types and momenta, from a sub-process passes the cuts...
Energy minKT(tcPDPtr p) const
Return the minimum allowed value of the transverse momentum of an outgoing parton.
double yHatMin() const
The minimum value of the rapidity of the hard sub-process (wrt.
Energy2 minTij(tcPDPtr pi, tcPDPtr po) const
Return the minimum allowed value of the negative of the squared invariant mass of an incoming parton ...
double theX2Min
The minimum value of the negative light-cone fraction of the hard sub-process.
Definition: Cuts.h:668
A concreate implementation of ClassDescriptionBase describing a concrete class with persistent data...
double theX1Max
The maximum value of the positive light-cone fraction of the hard sub-process.
Definition: Cuts.h:662
double theX1Min
The minimum value of the positive light-cone fraction of the hard sub-process.
Definition: Cuts.h:656
Energy2 maxScaleMin() const
Helper function used by the interface.
Ptr< JetFinder >::ptr theJetFinder
An optional jet finder used to define cuts on the level of reconstructed jets.
Definition: Cuts.h:710
Energy mHatMin() const
The minimum allowed value of .
Definition: Cuts.h:408
PersistentOStream is used to write objects persistently to a stream from which they can be read in ag...
double theYHatMin
The minimum value of the rapidity of the hard sub-process (wrt.
Definition: Cuts.h:644
vector< tcPDPtr > tcPDVector
A vector of transient pointers to const ParticleData objects.
Definition: Containers.h:42
void add(tOneCutPtr c)
Add a OneCutBase object.
Definition: Cuts.h:341
Energy2 maxS(const tcPDVector &pv) const
Return the maximum allowed value of the squared invariant mass of a set of outgoing partons of the gi...
TransientConstRCPtr is a simple wrapper around a bare const pointer which can be assigned to and from...
Definition: RCPtr.h:681
vector< typename Ptr< T >::transient_const_pointer > multiCutObjects() const
Return a vector of pointers to objects of the given class (with base class MultiCutBase).
double minYStar(tcPDPtr p) const
Return the minimum allowed rapidity of an outgoing parton of the given type in the center-of-mass sys...
TransientRCPtr is a simple wrapper around a bare pointer which can be assigned to and from an RCPtr a...
Definition: RCPtr.h:509
double theCurrentYHat
The total rapidity of hard sub-process (wrt.
Definition: Cuts.h:628
virtual bool initSubProcess(Energy2 shat, double yhat, bool mirror=false) const
Set information about the invariant mass squared, shat, and rapidity, yhat, of the hard sub-process...
double x1Max() const
The maximum value of the positive light-cone fraction of the hard sub-process.
const TwoCutVector & twoCuts() const
Return the objects defining cuts on pairs of particles in the hard sub-process.
Definition: Cuts.h:325
This is the main namespace within which all identifiers in ThePEG are declared.
Definition: FactoryBase.h:28
virtual void doinitrun()
Initialize this object.
Energy minKTClus(tcPDPtr pi, tcPDPtr pj) const
Return the minimum allowed value of the longitudinally invariant -algorithms distance measure...
Energy2 sHatMin() const
The minimum allowed value of .
Definition: Cuts.h:391
Energy2 theScaleMax
The maximum allowed value of the scale to be used in PDF&#39;s and coupling constants.
Definition: Cuts.h:686
Energy theMHatMin
The minimum allowed value of .
Definition: Cuts.h:633
bool theSubMirror
Set to true if a matrix element is should be using this cut and is mirrored along the z-axis ...
Definition: Cuts.h:716
Energy2 theSMax
The maximum allowed total invariant mass squared allowed for events to be considered.
Definition: Cuts.h:609
Cuts & operator=(const Cuts &)
The assignment operator is private and must never be called.
A SubProcess object represents a hard sub-process in a collision.
Definition: SubProcess.h:33
bool x2(double x) const
Check if the given is within the cuts.
QTY< 0, 1, 0 >::Type Energy
Energy.
Definition: Unitsystem.h:34
static ClassDescription< Cuts > initCuts
The static object used to initialize the description of this class.
Definition: Cuts.h:724
bool subMirror() const
Set true if a matrix element is should be using this cut and is mirrored along the z-axis ...
Definition: Cuts.h:487
double minRapidityMax(tcPDPtr p) const
Return the minimum allowed rapidity of an outgoing parton of the given type.
virtual IBPtr clone() const
Make a simple clone of this object.
Energy2 scaleMin() const
The minimum allowed value of the scale to be used in PDF&#39;s and coupling constants.
Definition: Cuts.h:470
vector< OneCutPtr > OneCutVector
A vector of OneCutBase pointers.
Definition: Cuts.h:59
bool sHat(Energy2 sh) const
Check if the given is within the cuts.
Definition: Cuts.h:401
virtual void describe() const
Describe the currently active cuts in the log file.
vector< tcPPtr > tcPVector
A vector of transient pointers to const Particle objects.
Definition: Containers.h:85
Energy2 scaleMax() const
The maximum allowed value of the scale to be used in PDF&#39;s and coupling constants.
Definition: Cuts.h:476
vector< typename Ptr< T >::transient_const_pointer > twoCutObjects() const
Return a vector of pointers to objects of the given class (with base class TwoCutBase).
const OneCutVector & oneCuts() const
Return the objects defining cuts on single outgoing partons from the hard sub-process.
Definition: Cuts.h:319
const MultiCutVector & multiCuts() const
Return the objects defining cuts on sets of outgoing particles from the hard sub-process.
Definition: Cuts.h:331
double minDeltaR(tcPDPtr pi, tcPDPtr pj) const
Return the minimum allowed value of of two outgoing partons of type pi and pj.
double yHatMax() const
The maximum value of the rapidity of the hard sub-process (wrt.
Cuts(Energy MhatMin=2 *GeV)
The default constructor.
MultiCutVector theMultiCuts
The objects defining cuts on sets of outgoing particles from the hard sub-process.
Definition: Cuts.h:704
Energy2 SMax() const
The maximum allowed total invariant mass squared allowed for events to be considered.
Definition: Cuts.h:362
double x2Max() const
The maximum value of the negative light-cone fraction of the hard sub-process.
Energy2 minScaleMax() const
Helper function used by the interface.
Energy mHatMax() const
The maximum allowed value of .
Definition: Cuts.h:413
double minX1Max() const
Helper function used by the interface.
ThePEG::Ptr< ParticleData >::transient_const_pointer tcPDPtr
Alias for a transient pointer to a const ParticleData .
Definition: Pointers.h:64
void add(tMultiCutPtr c)
Add a MultiCutBase object.
Definition: Cuts.h:351
double maxX1Min() const
Helper function used by the interface.
double minX2Max() const
Helper function used by the interface.
double x1Min() const
The minimum value of the positive light-cone fraction of the hard sub-process.
RCPtr is a reference counted (smart) pointer.
Definition: RCPtr.h:60
double minEta(tcPDPtr p) const
Return the minimum allowed pseudo-rapidity of an outgoing parton of the given type.
TwoCutVector theTwoCuts
The objects defining cuts on pairs of particles in the hard sub-process.
Definition: Cuts.h:698
Energy minMHatMax() const
Helper function used by the interface.
static void Init()
The standard Init function used to initialize the interfaces.
This is the decalaration of the Collision class.
Definition: Collision.h:34
The Interfaced class is derived from the InterfacedBase class adding a couple of things particular to...
Definition: Interfaced.h:38
Energy theMHatMax
The maximum allowed value of .
Definition: Cuts.h:638
vector< TwoCutPtr > TwoCutVector
A vector of TwoCutBase pointers.
Definition: Cuts.h:64
bool x1(double x) const
Check if the given is within the cuts.
double theYHatMax
The maximum value of the rapidity of the hard sub-process (wrt.
Definition: Cuts.h:650
double theX2Max
The maximum value of the negative light-cone fraction of the hard sub-process.
Definition: Cuts.h:674
bool scale(Energy2 Q2) const
Check if the given scale is within the cuts.
Definition: Cuts.h:481
void persistentOutput(PersistentOStream &os) const
Function used to write out object persistently.
The default concrete implementation of ClassTraitsBase.
Definition: ClassTraits.h:134
double maxX2Min() const
Helper function used by the interface.
const double epsilon
The smallest non-zero double.
Definition: Constants.h:63
Energy2 theCurrentSHat
The invariant mass squared of the hard sub-process of the event being considered. ...
Definition: Cuts.h:621
OneCutVector theOneCuts
The objects defining cuts on single outgoing partons from the hard sub-process.
Definition: Cuts.h:692
virtual ~Cuts()
The destructor.
Energy2 theScaleMin
The minimum allowed value of the scale to be used in PDF&#39;s and coupling constants.
Definition: Cuts.h:680
virtual void initialize(Energy2 smax, double Y)
Initialize this object specifying the maximum total invariant mass squared, smax, and the total rapid...
double minDurham(tcPDPtr pi, tcPDPtr pj) const
Return the minimum allowed value of the Durham -algorithms distance measure.
double maxYHatMin() const
Helper function used by the interface.
Ptr< JetFinder >::tptr jetFinder() const
Return the jet finder.
Definition: Cuts.h:336
BaseClassTraits describes the base classes of the templated class.
Definition: ClassTraits.h:156
double maxRapidityMin(tcPDPtr p) const
Return the maximum allowed rapidity of an outgoing parton of the given type.
vector< MultiCutPtr > MultiCutVector
A vector of MultiCutBase pointers.
Definition: Cuts.h:69
virtual IBPtr fullclone() const
Make a clone of this object, possibly modifying the cloned object to make it sane.
Energy2 minS(const tcPDVector &pv) const
Return the minimum allowed value of the squared invariant mass of a set of outgoing partons of the gi...
double x2Min() const
The minimum value of the negative light-cone fraction of the hard sub-process.
Energy2 sHatMax() const
The maximum allowed value of .
Definition: Cuts.h:396
vector< typename Ptr< T >::transient_const_pointer > oneCutObjects() const
Return a vector of pointers to objects of the given class (with base class OneCutBase).
Cuts is a class for implementing kinematical cuts in ThePEG.
Definition: Cuts.h:52
void add(tTwoCutPtr c)
Add a TwoCutBase object.
Definition: Cuts.h:346
The templated ClassTraitsBase class defines a set of default information about classes used by ThePEG...
Definition: ClassTraits.h:52
double maxYStar(tcPDPtr p) const
Return the minimum allowed rapidity of an outgoing parton of the given type in the center-of-mass sys...
virtual void initEvent()
Initialize this object for a new event.