ThePEG/PDF/PartonExtractor.h

0001 // -*- C++ -*-
0002 //
0003 // PartonExtractor.h is a part of ThePEG - Toolkit for HEP Event Generation
0004 // Copyright (C) 1999-2019 Leif Lonnblad
0005 //
0006 // ThePEG is licenced under version 3 of the GPL, see COPYING for details.
0007 // Please respect the MCnet academic guidelines, see GUIDELINES for details.
0008 //
0009 #ifndef ThePEG_PartonExtractor_H
0010 #define ThePEG_PartonExtractor_H
0011 // This is the declaration of the PartonExtractor class.
0012
0013 #include "ThePEG/Handlers/HandlerBase.h"
0014 #include "ThePEG/Handlers/LastXCombInfo.h"
0015 #include "ThePEG/PDF/PartonBin.h"
0016 #include "ThePEG/PDF/PartonBinInstance.h"
0017 #include "ThePEG/PDF/PDFBase.h"
0018 #include "ThePEG/PDT/ParticleData.h"
0019 #include "PartonExtractor.xh"
0020
0021 namespace ThePEG {
0022
0023 /**
0024  * The PartonExtractor is a base class defining the interface to
0025  * objects responsible for extracting partons from particles. It is
0026  * used by a SubProcessHandler which combines one PartonExtractor with
0027  * a number of MEBase objects which are the used in an XComb in a
0028  * StandardEventHandler to generate hard sub-processes.
0029  *
0030  * PartonExtractor inherits from the general HandlerBase class and
0031  * from the LastXCombInfo class to have easy acces to information
0032  * about the currently chosen hard sub-process.
0033  *
0034  * @see \ref PartonExtractorInterfaces "The interfaces"
0035  * defined for PartonExtractor.
0036  * @see SubProcessHandler
0037  * @see MEBase
0038  * @see EventHandler
0039  * @see StandardEventHandler
0040  * @see XComb
0041  * @see HandlerBase
0042  *
0043  */
0044 class PartonExtractor: public HandlerBase, public LastXCombInfo<> {
0045
0046   /** XComb is a friend. */
0047   friend class XComb;
0048
0049 public:
0050
0051   /** A map of PartonBinInstance objects indexed by the extracted parton. */
0052   typedef map<cPPtr,PBIPtr> PartonBinInstanceMap;
0053
0054 public:
0055
0056   /** @name Standard constructors and destructors. */
0057   //@{
0058   /**
0059    * Default constructor.
0060    */
0061   PartonExtractor();
0062
0063   /**
0064    * Destructor.
0065    */
0066   virtual ~PartonExtractor();
0067   //@}
0068
0069 public:
0070
0071   /** @name Virtual functions which may be overridden in sub-classes. */
0072   //@{
0073   /**
0074    * Return true if this parton extractor can handle the given types
0075    * of incoming particles.
0076    */
0077   virtual bool canHandle(const cPDPair &) { return true; }
0078
0079   /**
0080    * Return a vector of possible pairs of parton bins which can be
0081    * produced within a given maximum total particle-particle
0082    * invariant mass squared, \a maxEnergy sBin.
0083    */
0084   virtual PartonPairVec getPartons(Energy maxEnergy, const cPDPair &,
0085                    const Cuts &) const;
0086
0087   /**
0088    * May be overriden by sub-classes which have their own oppinion
0089    * about what scale to use in a hard subprocess. The default version
0090    * simply returns the previously selected scale.
0091    */
0092   virtual Energy2 newScale();
0093
0094   /**
0095    * Connect the remnants with the colour lines of the extracted
0096    * parton.
0097    */
0098   virtual void colourConnect(tPPtr particle, tPPtr parton,
0099                  const tPVector & remnants) const;
0100
0101   /**
0102    * If remnants has already been created for the given parton, remove
0103    * them from the given step and generate new remnants corresponding
0104    * to the parton newp and add them to the step. The new parton bins
0105    * are returned.
0106    * @throws Veto if remnant generation failed for whatever reason.
0107    */
0108   virtual PBIPair newRemnants(tPPair oldp, tPPair newp, tStepPtr step);
0109
0110   /**
0111    * Determine the number of random numbers needed to calculate
0112    * \f$\hat{s}\f$ and the product of all densitiy functions.
0113    */
0114   virtual pair<int,int> nDims(const PBPair & pbins);
0115
0116   /**
0117    * Prepare the given parton bin instances for generating a new
0118    * event.
0119    */
0120   virtual void prepare(const PBIPair & pbins);
0121
0122   /**
0123    * Update information on the given parton bin instances
0124    */
0125   virtual void updatePartonBinInstances(const PBIPair & pbins);
0126
0127   /**
0128    * Generate \f$l=\log(1/x)\f$ for all parton extractions.
0129    */
0130   virtual bool generateL(const PBIPair & pbins,
0131              const double * r1, const double * r2);
0132
0133   /**
0134    * Used by generateL() for each of the final parton
0135    * bins. Direction<0> is set to positive(negative) for the
0136    * first(second) final bin.
0137    */
0138   virtual void generateL(PartonBinInstance & pb, const double * r);
0139
0140   /**
0141    * Generate the rest of the degrees of freedom to calculate
0142    * \f$\hat{s}\f$ and the product of all densitiy functions.
0143    */
0144   virtual Energy2 generateSHat(Energy2 s, const PBIPair & pbins,
0145                    const double * r1, const double * r2,
0146                    bool mepartons = false);
0147
0148   /**
0149    * Return the product of all density functions. If noLastPDF.first
0150    * (.second) is true, then the PDF value multiplied by the momentum
0151    * fraction for the last extracted parton from the first (second)
0152    * incoming particle will be excluded.
0153    */
0154   virtual double fullFn(const PBIPair & pbins, Energy2 scale,
0155             pair<bool,bool> noLastPDF = make_pair(false,false));
0156
0157   /**
0158    * Construct remnants and add them to the step.
0159    */
0160   virtual void construct(const PBIPair & pbins, tStepPtr step) const;
0161
0162   /**
0163    * Construct remnants for partons created outside of this
0164    * extractor. Information about the incoming partons should be set
0165    * in \a pbins and the hard subprocess should be present in \a
0166    * sub. Generated remnants will be added to the \a step.
0167    * @throws Veto if remnant generation failed for whatever reason.
0168    */
0169   virtual void constructRemnants(const PBIPair & pbins, tSubProPtr sub,
0170              tStepPtr step) const;
0171
0172   /**
0173    * Get boost for hard subsystem and boost remnants. To be called
0174    * after re-constructing remnants and obtaining new momenta of the
0175    * partons entering the hard subsystem, but ignoring detailed energy
0176    * and momentum conservation. Perform boosts of the remnants to
0177    * conserve energy and momentum and return the boost needed for the
0178    * hard subsystem. \a bins contains the current state of the
0179    * incoming partons, including the momenta obtained after the
0180    * remnant generation. \a k1 and \a k2 contains the momenta of the
0181    * incoming partons before the remnant generation. If either side
0182    * has no new remnants, \a side1 and/or \a side2 should be false.
0183    */
0184   virtual LorentzRotation
0185   boostRemnants(PBIPair & bins, LorentzMomentum k1, LorentzMomentum k2,
0186        bool side1, bool side2) const;
0187   //@}
0188
0189   /** @name Access information about the current paron extraction. */
0190   //@{
0191   /**
0192    * Return the corresponding parton bin instance for a given
0193    * extracted parton.
0194    */
0195   tPBIPtr partonBinInstance(tcPPtr) const;
0196
0197   /**
0198    * Set the XComb object describing the current hard sub-process.
0199    */
0200   void select(tXCombPtr newXComb);
0201
0202   //@}
0203
0204   /**
0205    * The maximum number of attempts allowed when trying to generate
0206    * remnants.
0207    */
0208   int maxTries() const { return theMaxTries; }
0209
0210   /**
0211    * Return the PDFBase object to be used for the incoming particle
0212    * type. If one of theSpecialDensities matches the particle type it
0213    * is returned, otherwise if particle is a BeamParticleData use the
0214    * PDFBase object specified there. If also this fails, return a
0215    * NoPDF object.
0216    */
0217   tcPDFPtr getPDF(tcPDPtr particle) const;
0218
0219 protected:
0220
0221   /** @name Functions used by the main virtual functions. Some of
0222       these may be overridden in sub-classes. */
0223   //@{
0224
0225   /**
0226    * Used by generateSHat() for each of the final parton
0227    * bins. Direction<0> is set to positive(negative) for the
0228    * first(second) final bin, \a pb. Should ask the remnant handler to
0229    * generate what is needed to construct the extracted parton
0230    * momentum. \a r is a pointer to an array of random numbers to be
0231    * used and \a shat is the approximate invariant mass squared of the
0232    * hard system produced by the extracted parton and the primary
0233    * parton from the other side. \a first is the momentum of the
0234    * original incoming particle.
0235    *
0236    * @return false if no remnants could be generated.
0237    */
0238   virtual bool generate(PartonBinInstance & pb, const double * r,
0239             Energy2 shat, const Lorentz5Momentum & first,
0240             bool haveMEPartons = false);
0241
0242   /**
0243    * Used by the public fullFn() for each of the final parton bins.
0244    */
0245   virtual double fullFn(const PartonBinInstance & pb,
0246             bool noLastPDF = false);
0247
0248   /**
0249    * Used by the public construct() for each of the final parton
0250    * bins. If boost is false, no boost is necessary to give the
0251    * remnants proper momenta.
0252    */
0253   virtual void construct(PartonBinInstance & pb,
0254              tStepPtr step, bool boost = true) const;
0255
0256   /**
0257    * Used by the public newRemnants() for each of the parton bins.
0258    * @throws Veto if remnant generation failed for whatever reason.
0259    */
0260   PBIPtr newRemnants(tPBIPtr oldpb, tPPtr newp, const LorentzMomentum & k);
0261
0262   /**
0263    * Used by the public newRemnants() for each of the parton bins.
0264    */
0265   void addNewRemnants(tPBIPtr oldpb, tPBIPtr newpb, tStepPtr step);
0266
0267   /**
0268    * Transform remnant momentum. Assuming remnants have been generated
0269    * with momentum \a Pr without considering energy-momentum
0270    * conservation, shift the momentum, possibly compensating with the
0271    * momentum of the hard subsystem, \a Ph. For information the
0272    * momentum of the parton entering the hard subsystem from the other
0273    * side, \a k, and the momentum of the remnants parent particle , \a
0274    * P is given. Note that Direction<0> must be set to determine if
0275    * the parent particle is to be assumed to go in the positive or
0276    * negative direction.
0277    */
0278   virtual void transformRemnants(LorentzMomentum & Ph, LorentzMomentum & Pr,
0279                  const LorentzMomentum & k,
0280                  const LorentzMomentum & P) const;
0281
0282   /**
0283    * Construct remnants recursively for the parton represented by \a
0284    * pb. Used by constructRemnants(const PBIPair &, tSubProPtr, tStepPtr).
0285    * Shift the momentum, \a Ph, of the hard subsystem to conserve
0286    * energy and momentum if necessary. The momentum, \a k, of the
0287    * parton coming into the hard subsystem from the other side is
0288    * given for information. Note that Direction<0> must be set to
0289    * determine if the parent particle is to be assumed to go in the
0290    * positive or negative direction.
0291    * @throws Veto if remnant generation failed for whatever reason.
0292    */
0293   virtual void
0294   constructRemnants(PartonBinInstance & pb, LorentzMomentum & Ph,
0295             const LorentzMomentum & k) const;
0296   //@}
0297
0298 public:
0299
0300   /** @name Functions used by the persistent I/O system. */
0301   //@{
0302   /**
0303    * Function used to write out object persistently.
0304    * @param os the persistent output stream written to.
0305    */
0306   void persistentOutput(PersistentOStream & os) const;
0307
0308   /**
0309    * Function used to read in object persistently.
0310    * @param is the persistent input stream read from.
0311    * @param version the version number of the object when written.
0312    */
0313   void persistentInput(PersistentIStream & is, int version);
0314   //@}
0315
0316   /**
0317    * Standard Init function used to initialize the interface.
0318    */
0319   static void Init();
0320
0321 protected:
0322
0323   /**
0324    * Add parton bins to pbins for the given incoming particle and the
0325    * specified cuts.
0326    */
0327   virtual void addPartons(tPBPtr incoming ,const PDFCuts & cuts,
0328               tcPDFPtr pdf ,PartonVector & pbins) const;
0329
0330   /**
0331    * The NoPDF object.
0332    */
0333   tcPDFPtr noPDF() const { return theNoPDF; }
0334
0335   /**
0336    * Connect the first (\a anti) coloured particle in the given range
0337    * (not equal to \a parton) and connect it to the colour \a line.
0338    */
0339   template <typename Iterator>
0340   void findConnect(tColinePtr line, tPPtr parton, bool anti,
0341            Iterator first, Iterator last) const {
0342     for ( ; first != last; ++first ) {
0343       if ( *first != parton &&  (**first).hasColour(anti) &&
0344        !(**first).colourLine(anti) ) {
0345     line->addColoured(*first, anti);
0346     return;
0347       }
0348     }
0349     throw RemColException(*this);
0350   }
0351
0352 protected:
0353
0354   /** @name Clone Methods. */
0355   //@{
0356   /**
0357    * Make a simple clone of this object.
0358    * @return a pointer to the new object.
0359    */
0360   virtual IBPtr clone() const;
0361
0362   /** Make a clone of this object, possibly modifying the cloned object
0363    * to make it sane.
0364    * @return a pointer to the new object.
0365    */
0366   virtual IBPtr fullclone() const;
0367   //@}
0368
0369   /** @name Standard Interfaced functions. */
0370   //@{
0371
0372   /**
0373    * Finalize this object. Called in the run phase just after a
0374    * run has ended. Used eg. to write out statistics.
0375    */
0376   virtual void dofinish();
0377   //@}
0378
0379 private:
0380
0381   /**
0382    * The PartonBinInstance's used mapped to the respective partons.
0383    */
0384   PartonBinInstanceMap& partonBinInstances() const {
0385     assert(lastXCombPtr());
0386     return lastXCombPtr()->partonBinInstanceMap();
0387   }
0388
0389   /**
0390    * A list of special PDFBase objects to be used.
0391    */
0392   vector<PDFPtr> theSpecialDensities;
0393
0394   /**
0395    *  PDFBase object to override  first PDF
0396    */
0397   PDFPtr theFirstPDF;
0398
0399   /**
0400    *  PDFBase object to override second PDF
0401    */
0402   PDFPtr theSecondPDF;
0403
0404   /**
0405    * The NoPDF object.
0406    */
0407   PDFPtr theNoPDF;
0408
0409   /**
0410    * The maximum number of tries allowed when trying to produce
0411    * remnants.
0412    */
0413   int theMaxTries;
0414
0415   /**
0416    * True if this extractor should override the \f$l\f$-generation in
0417    * the PDFs and generate a flat distribution in \f$\log(\hat{s})\f$
0418    * and y.
0419    */
0420   bool flatSHatY;
0421
0422 private:
0423
0424   /**
0425    * Describe a concrete class with persistent data.
0426    */
0427   static ClassDescription<PartonExtractor> initPartonExtractor;
0428
0429   /**
0430    *  Private and non-existent assignment operator.
0431    */
0432   PartonExtractor & operator=(const PartonExtractor &) = delete;
0433
0434 };
0435
0436 /** @cond TRAITSPECIALIZATIONS */
0437
0438 /** This template specialization informs ThePEG about the base classes
0439  *  of PartonExtractor. */
0440 template <>
0441 struct BaseClassTrait<PartonExtractor,1>: public ClassTraitsType {
0442   /** Typedef of the first base class of PartonExtractor. */
0443   typedef HandlerBase NthBase;
0444 };
0445
0446 /** This template specialization informs ThePEG about the name of the
0447  *  PartonExtractor class. */
0448 template <>
0449   /** Return a platform-independent class name */
0450 struct ClassTraits<PartonExtractor>: public ClassTraitsBase<PartonExtractor> {
0451   static string className() { return "ThePEG::PartonExtractor"; }
0452 };
0453
0454 /** @endcond */
0455
0456 }
0457
0458 #endif /* ThePEG_PartonExtractor_H */