PerturbativeDecayer.h

// -*- C++ -*-
#ifndef Herwig_PerturbativeDecayer_H
#define Herwig_PerturbativeDecayer_H
//
// This is the declaration of the PerturbativeDecayer class.
//
 
#include "Herwig/Decay/DecayIntegrator.h"
#include "Herwig/Shower/ShowerAlpha.h"
#include "Herwig/Shower/ShowerInteraction.h"
 
namespace Herwig {
 
using namespace ThePEG;
 
class PerturbativeDecayer: public DecayIntegrator {
 
protected:
  
  enum dipoleType {FFa, FFc, IFa, IFc, IFba, IFbc, FFg};
 
  enum phaseSpaceRegion {emissionFromB,emissionFromC,emissionFromA1,emissionFromA2,deadZone};
 
  struct DipoleType {
 
    DipoleType() {}
 
    DipoleType(dipoleType a, ShowerInteraction b)
      : type(a), interaction(b)
    {}
    
    dipoleType type;
 
    ShowerInteraction interaction;
  };
 
public:
 
  PerturbativeDecayer() : inter_(ShowerInteraction::QCD),
                          pTmin_(GeV), useMEforT2_(true),
                          C_(5.), ymax_(10.), phaseOpt_(1),
                          pT_(ZERO),mb_(ZERO), e_(0.),
                          s_(0.), e2_(0.), s2_(0.), enhance_(1.)
  {}
 
  virtual POWHEGType hasPOWHEGCorrection() {return No;}
 
  virtual RealEmissionProcessPtr generateHardest(RealEmissionProcessPtr);
 
  virtual RealEmissionProcessPtr applyHardMatrixElementCorrection(RealEmissionProcessPtr);
 
public:
 
  void persistentOutput(PersistentOStream & os) const;
 
  void persistentInput(PersistentIStream & is, int version);
 
  static void Init();
 
protected:
 
  virtual double matrixElementRatio(const Particle & inpart, const ParticleVector & decay2,
                                    const ParticleVector & decay3, MEOption meopt,
                                    ShowerInteraction inter);
 
  bool identifyDipoles(vector<DipoleType> & dipoles,
                       PPtr & aProgenitor,
                       PPtr & bProgenitor,
                       PPtr & cProgenitor,
                       ShowerInteraction inter) const;
 
  ShowerAlphaPtr alphaS() {return alphaS_;}
 
  ShowerAlphaPtr alphaEM() {return alphaEM_;}
 
  vector<Lorentz5Momentum> hardMomenta(PPtr in, PPtr emitter, 
                                       PPtr spectator, 
                                       const vector<DipoleType>  & dipoles,
                                       int i, bool inDeadZone);
 
  bool calcMomenta(int j, Energy pT, double y, double phi, double& xg, 
                   double& xs, double& xe, double& xe_z, 
                   vector<Lorentz5Momentum>& particleMomenta);
 
  bool psCheck(const double xg, const double xs);
 
  pair<double,double> calculateDipole(const DipoleType & dipoleId,
                                      const Particle & inpart,
                                      const ParticleVector & decay3);
 
  double dipoleSpinFactor(tcPDPtr emitter, double z);
 
  double colourCoeff(tcPDPtr emitter, tcPDPtr spectator,
                     tcPDPtr other, DipoleType dipole);
  
  void getColourLines(RealEmissionProcessPtr real);
 
  RealEmissionProcessPtr getHardEvent(RealEmissionProcessPtr born,
                                      bool inDeadZone,
                                      ShowerInteraction inter);
 
  bool inTotalDeadZone(double xg, double xs,
                       const vector<DipoleType>  & dipoles,
                       int i);
 
  phaseSpaceRegion inInitialFinalDeadZone(double xg, double xa, double a, double c) const;
 
  phaseSpaceRegion inFinalFinalDeadZone(double xb, double xc, double b, double c) const;
 
  bool useMEforT2() const {return useMEforT2_;}
 
protected:
 
  const Energy & pT() const { return pT_;}
 
  const Energy & mb() const {return mb_;}
 
  const double & e() const {return e_;}
 
  const double & s() const {return s_;}
 
  const double & e2() const {return e2_;}
 
  const double & s2() const {return s2_;}
 
private:
 
  PerturbativeDecayer & operator=(const PerturbativeDecayer &) = delete;
 
private:
 
  ShowerInteraction inter_;
  
  ShowerAlphaPtr alphaS_;
  
  ShowerAlphaPtr alphaEM_;
 
  Energy pTmin_;
 
  bool useMEforT2_;
 
  double C_;
 
  double ymax_;
 
  unsigned int phaseOpt_;
 
private:
 
  Energy pT_;
 
  Energy mb_;
 
  double e_;
 
  double s_;
 
  double e2_;
 
  double s2_;
 
  mutable double enhance_;
};
 
}
 
#endif /* Herwig_PerturbativeDecayer_H */