Page Menu
Home
HEPForge
Search
Configure Global Search
Log In
Files
F8725591
No One
Temporary
Actions
View File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Flag For Later
Size
22 KB
Subscribers
None
View Options
diff --git a/Handlers/StdDependentXComb.cc b/Handlers/StdDependentXComb.cc
--- a/Handlers/StdDependentXComb.cc
+++ b/Handlers/StdDependentXComb.cc
@@ -1,165 +1,167 @@
// -*- C++ -*-
//
// StdDependentXComb.cc is a part of ThePEG - Toolkit for HEP Event Generation
// Copyright (C) 1999-2007 Leif Lonnblad
// Copyright (C) 2009-2010 Simon Platzer
//
// ThePEG is licenced under version 2 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
//
// This is the implementation of the non-inlined, non-templated member
// functions of the StdDependentXComb class.
//
#include "StdDependentXComb.h"
#include "StdXCombGroup.h"
#include "ThePEG/MatrixElement/MEGroup.h"
#include "ThePEG/Cuts/Cuts.h"
#include "ThePEG/Handlers/CascadeHandler.h"
#include "ThePEG/Persistency/PersistentOStream.h"
#include "ThePEG/Persistency/PersistentIStream.h"
using namespace ThePEG;
StdDependentXComb::StdDependentXComb(tStdXCombPtr newHead,
const PBPair & newPartonBins, tMEPtr newME,
const DiagramVector & newDiagrams)
: StandardXComb(newHead->maxEnergy(),newHead->particles(),
newHead->eventHandlerPtr(),
const_ptr_cast<tSubHdlPtr>(newHead->subProcessHandler()),
newHead->pExtractor(),newHead->CKKWHandler(),
newPartonBins,newHead->cuts(),newME,newDiagrams,newHead->mirror(),
newHead), resetIncoming(true) {}
StdDependentXComb::StdDependentXComb()
: StandardXComb(), resetIncoming(true) {}
StdDependentXComb::~StdDependentXComb() { }
tSubProPtr StdDependentXComb::construct() {
// first get the meMomenta in their CMS, as this may
// not be the case
Boost cm = (meMomenta()[0] + meMomenta()[1]).findBoostToCM();
if ( cm.mag2() > Constants::epsilon ) {
for ( vector<Lorentz5Momentum>::iterator m = meMomenta().begin();
m != meMomenta().end(); ++m ) {
*m = m->boost(cm);
}
}
tSubProPtr sub = StandardXComb::construct();
sub->head(head()->subProcess());
sub->groupWeight(lastCrossSection()/head()->lastCrossSection());
return sub;
}
void StdDependentXComb::setProcess() {
meMomenta().resize(mePartonData().size());
}
void StdDependentXComb::setPartonBinInstances(Energy2 scale) {
PBIPair newBins;
Direction<0> dir(true);
newBins.first =
new_ptr(PartonBinInstance(lastPartons().first,partonBins().first,scale));
dir.reverse();
newBins.second =
new_ptr(PartonBinInstance(lastPartons().second,partonBins().second,scale));
resetPartonBinInstances(newBins);
setPartonBinInfo();
lastPartons().first->scale(partonBinInstances().first->scale());
lastPartons().second->scale(partonBinInstances().second->scale());
}
void StdDependentXComb::setIncomingPartons() {
if ( !resetIncoming )
return;
resetIncoming = false;
clean();
createPartonBinInstances();
setPartonBinInfo();
lastParticles(head()->lastParticles());
lastPartons(make_pair(mePartonData()[0]->produceParticle(Lorentz5Momentum()),
mePartonData()[1]->produceParticle(Lorentz5Momentum())));
Lorentz5Momentum pFirst = meMomenta()[0];
Lorentz5Momentum pSecond = meMomenta()[1];
if ( head()->matrixElement()->wantCMS() ) {
Boost toLab = (head()->lastPartons().first->momentum() +
head()->lastPartons().second->momentum()).boostVector();
if ( toLab.mag2() > Constants::epsilon ) {
pFirst.boost(toLab);
pSecond.boost(toLab);
}
}
lastPartons().first->set5Momentum(pFirst);
lastPartons().second->set5Momentum(pSecond);
lastS((lastParticles().first->momentum() +
lastParticles().second->momentum()).m2());
lastSHat((lastPartons().first->momentum() +
lastPartons().second->momentum()).m2());
lastP1P2(make_pair(0.0, 0.0));
double x1 =
lastPartons().first->momentum().plus()/
lastParticles().first->momentum().plus();
double x2 =
lastPartons().second->momentum().minus()/
lastParticles().second->momentum().minus();
lastX1X2(make_pair(x1,x2));
lastY(log(lastX1()/lastX2())*0.5);
}
CrossSection StdDependentXComb::dSigDR() {
setProcess();
setIncomingPartons();
lastScale(head()->lastScale());
lastAlphaS(head()->lastAlphaS());
lastAlphaEM(head()->lastAlphaEM());
- if ( (!willPassCuts() && !matrixElement()->headCuts()) ||
+ if ( (!willPassCuts() &&
+ !matrixElement()->headCuts() &&
+ !matrixElement()->ignoreCuts()) ||
!matrixElement()->apply() ) {
subProcess(SubProPtr());
lastCrossSection(ZERO);
return ZERO;
}
lastPDFWeight(head()->lastPDFWeight());
CrossSection xsec = matrixElement()->dSigHatDR() * lastPDFWeight();
subProcess(SubProPtr());
lastCrossSection(xsec);
return xsec;
}
void StdDependentXComb::Init() {}
void StdDependentXComb::persistentOutput(PersistentOStream &) const {
}
void StdDependentXComb::persistentInput(PersistentIStream &, int) {
}
ClassDescription<StdDependentXComb> StdDependentXComb::initStdDependentXComb;
diff --git a/MatrixElement/MEBase.h b/MatrixElement/MEBase.h
--- a/MatrixElement/MEBase.h
+++ b/MatrixElement/MEBase.h
@@ -1,584 +1,590 @@
// -*- C++ -*-
//
// MEBase.h is a part of ThePEG - Toolkit for HEP Event Generation
// Copyright (C) 1999-2011 Leif Lonnblad
//
// ThePEG is licenced under version 2 of the GPL, see COPYING for details.
// Please respect the MCnet academic guidelines, see GUIDELINES for details.
//
#ifndef ThePEG_MEBase_H
#define ThePEG_MEBase_H
// This is the declaration of the MEBase class.
#include "ThePEG/Handlers/HandlerBase.h"
#include "ThePEG/EventRecord/SubProcess.h"
#include "ThePEG/MatrixElement/DiagramBase.h"
#include "ThePEG/MatrixElement/ColourLines.h"
#include "ThePEG/MatrixElement/Amplitude.h"
#include "ThePEG/Handlers/LastXCombInfo.h"
#include "ThePEG/Handlers/StandardXComb.fh"
#include "ReweightBase.h"
#include "MEBase.fh"
namespace ThePEG {
/**
* The MEBase class is the base class of all objects
* representing hard matrix elements in ThePEG. There are three
* methods which must be overridden by a concrete subclass:<BR>
*
* includedDiagrams(tcPDPair) should return a vector of DiagramBase
* objects describing the diagrams used for this matrix element for
* the given pair of incoming parton types. These DiagramBases are
* used to identify the incoming and outgoing partons which can be
* handled by the process generation scheme, and is also used to
* cnstruct a corresponding SubProcess object.
*
* scale() should return the scale associated with the phase space
* point set with the last call to setKinematics(...) or
* generateKinematics(...).
*
* me() should return the the matrix element squared using the the
* type and momentum of the incoming and outgoing partons, previously
* set by the setKinematics(...) or generateKinematics(...) member
* functions, accessible through the methods meMomenta() and
* mePartonData() inherited from LastXCombInfo, and/or from
* information stored by sub classes. The returned value should be
* dimensionless suitable scaled by the total invariant mass squared
* (accessible through the sHat() member function). Any user of this
* method must make sure that the setKinematics(...) member function
* has been appropriately called before.
*
* colourGeometries() should return a Selector with the possible
* ColourLines objects weighted by their relative probabilities given
* the information set by the last call to setKinematics(...) or
* generateKinematics(...).
*
* There are other virtula functions which may be overridden as listed
* below.
*
* @see \ref MEBaseInterfaces "The interfaces"
* defined for MEBase.
* @see DiagramBase
* @see ColourLines
*
*/
class MEBase: public HandlerBase, public LastXCombInfo<StandardXComb> {
public:
/** A vector of pointers to DiagramBase objects. */
typedef vector<DiagPtr> DiagramVector;
/** The size_type used in the DiagramVector. */
typedef DiagramVector::size_type DiagramIndex;
/** A vector of pointers to ReweightBase objects. */
typedef vector<ReweightPtr> ReweightVector;
public:
/** @name Standard constructors and destructors. */
//@{
/**
* Default constructor.
*/
MEBase();
/**
* Destructor.
*/
virtual ~MEBase();
//@}
public:
/** @name Virtual functions to be overridden by sub-classes.. */
//@{
/**
* Return the order in \f$\alpha_S\f$ in which this matrix element
* is given.
*/
virtual unsigned int orderInAlphaS() const = 0;
/**
* Return the order in \f$\alpha_{EM}\f$ in which this matrix
* element is given. Returns 0.
*/
virtual unsigned int orderInAlphaEW() const = 0;
/**
* Return the matrix element for the kinematical configuation
* previously provided by the last call to setKinematics(), suitably
* scaled by sHat() to give a dimension-less number.
*/
virtual double me2() const = 0;
/**
* Return the scale associated with the phase space point provided
* by the last call to setKinematics().
*/
virtual Energy2 scale() const = 0;
/**
* Return the value of \f$\alpha_S\f$ associated with the phase
* space point provided by the last call to setKinematics(). This
* versions returns SM().alphaS(scale()).
*/
virtual double alphaS() const;
/**
* Return the value of \f$\alpha_EM\f$ associated with the phase
* space point provided by the last call to setKinematics(). This
* versions returns SM().alphaEM(scale()).
*/
virtual double alphaEM() const;
/**
* Set the typed and momenta of the incoming and outgoing partons to
* be used in subsequent calls to me() and colourGeometries().
*/
void setKinematics(tPPair in, const PVector & out);
/**
* Set the typed and momenta of the incoming and outgoing partons to
* be used in subsequent calls to me() and colourGeometries()
* according to the associated XComb object. If the function is
* overridden in a sub class the new function must call the base
* class one first.
*/
virtual void setKinematics();
/**
* construct the spin information for the interaction
*/
virtual void constructVertex(tSubProPtr sub);
/**
* The number of internal degreed of freedom used in the matrix
* element. This default version returns 0;
*/
virtual int nDim() const;
/**
* Generate internal degrees of freedom given nDim() uniform random
* numbers in the interval ]0,1[. To help the phase space generator,
* the 'dSigHatDR' should be a smooth function of these numbers,
* although this is not strictly necessary. The return value should
* be true of the generation succeeded. If so the generated momenta
* should be stored in the meMomenta() vector.
*/
virtual bool generateKinematics(const double * r) = 0;
/**
* Return true, if this matrix element expects
* the incoming partons in their center-of-mass system
*/
virtual bool wantCMS() const { return true; }
/**
* If this is a dependent matrix element in a ME group, return true,
* if cuts should be inherited from the head matrix element, i.e. no
* cut is being applied to the dependent matrix element if the head
* configuration has passed the cuts.
*/
virtual bool headCuts() const { return false; }
/**
* If this is a dependent matrix element in a ME group, return true,
+ * if cuts should be ignored.
+ */
+ virtual bool ignoreCuts() const { return false; }
+
+ /**
+ * If this is a dependent matrix element in a ME group, return true,
* if it applies to the process set in lastXComb()
*/
virtual bool apply() const { return true; }
/**
* Return the matrix element squared differential in the variables
* given by the last call to generateKinematics().
*/
virtual CrossSection dSigHatDR() const = 0;
/**
* Return true, if this matrix element will generate momenta for the
* incoming partons itself. The matrix element is required to store
* the incoming parton momenta in meMomenta()[0,1]. No mapping in
* tau and y is performed by the PartonExtractor object, if a
* derived class returns true here. The phase space jacobian is to
* include a factor 1/(x1 x2).
*/
virtual bool haveX1X2() const { return false; }
/**
* Return true, if this matrix element provides the PDF
* weight for the first incoming parton itself.
*/
virtual bool havePDFWeight1() const { return false; }
/**
* Return true, if this matrix element provides the PDF
* weight for the second incoming parton itself.
*/
virtual bool havePDFWeight2() const { return false; }
/**
* Return true, if the XComb steering this matrix element
* should keep track of the random numbers used to generate
* the last phase space point
*/
virtual bool keepRandomNumbers() const { return false; }
/**
* Comlete a SubProcess object using the internal degrees of freedom
* generated in the last generateKinematics() (and possible other
* degrees of freedom which was intergated over in dSigHatDR(). This
* default version does nothing. Will be made purely virtual in the
* future.
*/
virtual void generateSubCollision(SubProcess &);
/**
* Clear the information previously provided by a call to
* setKinematics(...).
*/
virtual void clearKinematics();
/**
* Add all possible diagrams with the add() function.
*/
virtual void getDiagrams() const = 0;
/**
* Return true, if this matrix element does not want to
* make use of mirroring processes; in this case all
* possible partonic subprocesses with a fixed assignment
* of incoming particles need to be provided through the diagrams
* added with the add(...) method.
*/
virtual bool noMirror () const { return false; }
/**
* Return all possible diagrams.
*/
const DiagramVector & diagrams() const {
if ( theDiagrams.empty() ) getDiagrams();
return theDiagrams;
}
/**
* Return a Selector with possible colour geometries for the selected
* diagram weighted by their relative probabilities.
*/
virtual Selector<const ColourLines *>
colourGeometries(tcDiagPtr diag) const = 0;
/**
* Select a ColpurLines geometry. The default version returns a
* colour geometry selected among the ones returned from
* colourGeometries(tcDiagPtr).
*/
virtual const ColourLines &
selectColourGeometry(tcDiagPtr diag) const;
/**
* With the information previously supplied with the
* setKinematics(...) method, a derived class may optionally
* override this method to weight the given diagrams with their
* (although certainly not physical) relative probabilities.
*/
virtual Selector<DiagramIndex> diagrams(const DiagramVector &) const {
return Selector<DiagramIndex>();
}
/**
* Select a diagram. Default version uses diagrams(const
* DiagramVector &) to select a diagram according to the
* weights. This is the only method used that should be outside of
* MEBase.
*/
virtual DiagramIndex diagram(const DiagramVector &) const;
/**
* Return true if this matrix element has associated (p)reWeight
* objects assigned.
*/
inline bool reweighted() const {
return reweights.size() > 0 || preweights.size() > 0;
}
/**
* With the information previously supplied with the
* setKinematics(...) methods, return the combined effects of the
* reweighters.
*/
double reWeight() const;
/**
* With the information previously supplied with the
* setKinematics(...) methods, return the comined effects of the
* peweighters.
*/
double preWeight() const;
/**
* Add objects to the list of reweighters.
*/
void addReweighter(tReweightPtr rw);
/**
* Add objects to the list of preweighters.
*/
void addPreweighter(tReweightPtr rw);
/**
* Return the amplitude associated with this matrix element. This
* function is allowed to return the null pointer if the amplitude
* is not available.
*/
Ptr<Amplitude>::pointer amplitude() const { return theAmplitude; }
/**
* Set the amplitude associated with this matrix element.
*/
void amplitude(Ptr<Amplitude>::pointer amp) { theAmplitude = amp; }
//@}
public:
/** @name Acces information about the last generated phase space point. */
//@{
/**
* Return the last set invariant mass squared.
*/
Energy2 sHat() const { return theLastSHat; }
/**
* Return the factor with which this matrix element was last
* pre-weighted.
*/
double preweight() const { return lastPreweight; }
/**
* Inform this matrix element that a new phase space
* point is about to be generated, so all caches should
* be flushed.
*/
virtual void flushCaches() {}
/**
* Set the XComb object to be used in the next call to
* generateKinematics() and dSigHatDR().
*/
virtual void setXComb(tStdXCombPtr);
/**
* Retrieve information obtained in the calculation of the cross
* section to be used later when selecting diagrams and colour flow.
*/
const DVector & meInfo() const;
/**
* Save information obtained in the calculation of the cross
* section to be used later when selecting diagrams and colour flow.
*/
void meInfo(const DVector & info) const;
/**
* If this matrix element is to be used together with others for
* CKKW reweighting and veto, this should give the multiplicity of
* outgoing particles in the highest multiplicity matrix element in
* the group.
*/
int maxMultCKKW() const { return theMaxMultCKKW; }
/**
* If this matrix element is to be used together with others for
* CKKW reweighting and veto, this should give the multiplicity of
* outgoing particles in the lowest multiplicity matrix element in
* the group.
*/
int minMultCKKW() const { return theMinMultCKKW; }
/**
* Set veto scales on the particles at the given
* SubProcess which has been generated using this
* matrix element.
*/
virtual void setVetoScales(tSubProPtr) const {}
//@}
public:
/** @name Functions used by the persistent I/O system. */
//@{
/**
* Function used to write out object persistently.
* @param os the persistent output stream written to.
*/
void persistentOutput(PersistentOStream & os) const;
/**
* Function used to read in object persistently.
* @param is the persistent input stream read from.
* @param version the version number of the object when written.
*/
void persistentInput(PersistentIStream & is, int version);
//@}
/**
* Standard Init function used to initialize the interfaces.
*/
static void Init();
protected:
/**
* To be used by sub classes in the getDiagrams() method to add
* included diagrams.
*/
void add(DiagPtr dp) const { theDiagrams.push_back(dp); }
/**
* Access the momenta set by the last call to generateKinematics().
*/
vector<Lorentz5Momentum> & meMomenta();
using LastXCombInfo<StandardXComb>::meMomenta;
/**
* Set the matrix element squared as calculated
* for the last phase space point. This may optionally
* be used by a matrix element for caching.
*/
void lastME2(double v) const;
using LastXCombInfo<StandardXComb>::lastME2;
/**
* Set the partonic cross section as calculated
* for the last phase space point. This may optionally
* be used by a matrix element for caching.
*/
void lastMECrossSection(CrossSection v) const;
using LastXCombInfo<StandardXComb>::lastMECrossSection;
/**
* Set the PDF weight as calculated
* for the last phase space point, if the matrix
* element does supply PDF weights. This may optionally
* be used by a matrix element for caching.
*/
void lastMEPDFWeight(double v) const;
using LastXCombInfo<StandardXComb>::lastMEPDFWeight;
/**
* Set the last jacobian obtained when generating the kinematics for
* the call to dSigHatDR.
*/
void jacobian(double j);
using LastXCombInfo<StandardXComb>::jacobian;
/**
* Initialize all member variables from another
* MEBase object.
*
* @TODO remove?
*/
void use(tcMEPtr other);
/**
* Initialize the diagrams from another MEBase object.
*/
void useDiagrams(tcMEPtr other) const;
private:
/**
* The diagrams included for this matrix element.
*/
mutable DiagramVector theDiagrams;
/**
* The last set invariant mass squared.
*/
Energy2 theLastSHat;
/**
* The reweight objects modifying this matrix element.
*/
ReweightVector reweights;
/**
* The preweight objects modifying this matrix element.
*/
ReweightVector preweights;
/**
* The factor with which this matrix element was last pre-weighted.
*/
mutable double lastPreweight;
/**
* The amplitude associated with this matrix element.
*/
Ptr<Amplitude>::pointer theAmplitude;
/**
* If this matrix element is to be used together with others for
* CKKW reweighting and veto, this should give the multiplicity of
* outgoing particles in the highest multiplicity matrix element in
* the group.
*/
int theMaxMultCKKW;
/**
* If this matrix element is to be used together with others for
* CKKW reweighting and veto, this should give the multiplicity of
* outgoing particles in the lowest multiplicity matrix element in
* the group.
*/
int theMinMultCKKW;
private:
/**
* Describe an abstract base class with persistent data.
*/
static AbstractClassDescription<MEBase> initMEBase;
/**
* Private and non-existent assignment operator.
*/
MEBase & operator=(const MEBase &);
};
}
namespace ThePEG {
/** @cond TRAITSPECIALIZATIONS */
/**
* This template specialization informs ThePEG about the base class of
* MEBase.
*/
template <>
struct BaseClassTrait<MEBase,1>: public ClassTraitsType {
/** Typedef of the base class of MEBase. */
typedef HandlerBase NthBase;
};
/**
* This template specialization informs ThePEG about the name of the
* MEBase class.
*/
template <>
struct ClassTraits<MEBase>: public ClassTraitsBase<MEBase> {
/** Return the class name. */
static string className() { return "ThePEG::MEBase"; }
};
/** @endcond */
}
#include "ThePEG/Handlers/StandardXComb.h"
#endif /* ThePEG_MEBase_H */
File Metadata
Details
Attached
Mime Type
text/x-diff
Expires
Tue, Jan 21, 1:59 AM (1 d, 17 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
4243580
Default Alt Text
(22 KB)
Attached To
rTHEPEGHG thepeghg
Event Timeline
Log In to Comment