herwig is hosted by Hepforge, IPPP Durham
Herwig  7.2.1
MatchboxPhasespace.h
1 // -*- C++ -*-
2 //
3 // MatchboxPhasespace.h is a part of Herwig - A multi-purpose Monte Carlo event generator
4 // Copyright (C) 2002-2019 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_MatchboxPhasespace_H
10 #define HERWIG_MatchboxPhasespace_H
11 //
12 // This is the declaration of the MatchboxPhasespace class.
13 //
14 
15 #include "ThePEG/Handlers/StandardXComb.h"
16 #include "ThePEG/Handlers/HandlerBase.h"
17 #include "ThePEG/MatrixElement/Tree2toNDiagram.h"
18 #include "Herwig/MatrixElement/Matchbox/Utility/LastMatchboxXCombInfo.h"
19 #include "Herwig/MatrixElement/Matchbox/Utility/ProcessData.fh"
20 #include "Herwig/MatrixElement/Matchbox/MatchboxFactory.fh"
21 #include "Herwig/MatrixElement/Matchbox/Phasespace/PhasespaceCouplings.h"
22 
23 namespace Herwig {
24 
25 using namespace ThePEG;
26 
34 struct StreamingRnd {
35 
39  const double* numbers;
40 
44  size_t nRnd;
45 
50  : numbers(0), nRnd(0) {}
51 
55  explicit StreamingRnd(const double* newNumbers,
56  size_t n)
57  : numbers(newNumbers), nRnd(n) {}
58 
62  inline double operator()() {
63  assert(numbers && nRnd > 0);
64  const double ret = numbers[0];
65  ++numbers; --nRnd;
66  return ret;
67  }
68 
69 };
70 
80  public HandlerBase,
81  public LastXCombInfo<StandardXComb>,
82  public LastMatchboxXCombInfo {
83 
84 public:
85 
92 
96  virtual ~MatchboxPhasespace();
98 
99 public:
100 
105  virtual void setXComb(tStdXCombPtr xc) {
106  theLastXComb = xc;
107  lastMatchboxXComb(xc);
108  }
109 
113  Ptr<MatchboxFactory>::tcptr factory() const;
114 
118  Ptr<ProcessData>::tptr processData() const;
119 
123  virtual double generateKinematics(const double* r,
124  vector<Lorentz5Momentum>& momenta);
125 
129  virtual double generateTwoToNKinematics(const double*,
130  vector<Lorentz5Momentum>& momenta) = 0;
131 
135  virtual double generateTwoToOneKinematics(const double*,
136  vector<Lorentz5Momentum>& momenta);
137 
142  virtual int nDim(const cPDVector&) const;
143 
148  virtual int nDimPhasespace(int nFinal) const = 0;
149 
154  virtual bool haveX1X2() const { return false; }
155 
160  virtual bool wantCMS() const { return true; }
161 
165  bool useMassGenerators() const { return theUseMassGenerators; }
166 
170  virtual Selector<MEBase::DiagramIndex> selectDiagrams(const MEBase::DiagramVector&) const;
171 
176  pair<double,Lorentz5Momentum> timeLikeWeight(const Tree2toNDiagram& diag,
177  int branch, double flatCut) const;
178 
183  double spaceLikeWeight(const Tree2toNDiagram& diag,
184  const Lorentz5Momentum& incoming,
185  int branch, double flatCut) const;
186 
190  double diagramWeight(const Tree2toNDiagram& diag) const {
191  assert( !diagramWeights().empty() );
192  return diagramWeights().find(diag.id())->second;
193  }
194 
198  void fillDiagramWeights(double flatCut = 0.0);
199 
204  diagramWeights().clear();
205  }
206 
210  Ptr<MatchboxPhasespace>::ptr cloneMe() const {
211  return dynamic_ptr_cast<Ptr<MatchboxPhasespace>::ptr>(clone());
212  }
213 
217  virtual void cloneDependencies(const std::string& prefix = "");
218 
219 public:
220 
224  virtual bool isInvertible() const { return false; }
225 
230  virtual double invertKinematics(const vector<Lorentz5Momentum>& momenta,
231  double* r) const;
232 
237  virtual double invertTwoToNKinematics(const vector<Lorentz5Momentum>&,
238  double*) const {
239  return 0.;
240  }
241 
246  virtual double invertTwoToOneKinematics(const vector<Lorentz5Momentum>&, double*) const;
247 
248 public:
249 
253  void singularLimit(size_t i, size_t j) {
254  if ( i > j )
255  swap(i,j);
256  singularLimits().insert(make_pair(i,j));
257  }
258 
262  const pair<size_t,size_t>& lastSingularIndices() const {
263  assert(lastSingularLimit() != singularLimits().end());
264  return *lastSingularLimit();
265  }
266 
270  bool matchConstraints(const vector<Lorentz5Momentum>& momenta);
271 
272 protected:
273 
279  void setCoupling(long a, long b, long c,
280  double coupling, bool includeCrossings = true);
281 
282 public:
283 
290  void persistentOutput(PersistentOStream & os) const;
291 
297  void persistentInput(PersistentIStream & is, int version);
299 
300 public:
301 
308  static void Init();
309 
310 
311 // If needed, insert declarations of virtual function defined in the
312 // InterfacedBase class here (using ThePEG-interfaced-decl in Emacs).
313 
314 
315 private:
316 
321 
326 
330  Ptr<PhasespaceCouplings>::ptr theCouplings;
331 
335  string doSetCoupling(string);
336 
340  string doSetPhysicalCoupling(string);
341 
348 
355 
360  MatchboxPhasespace & operator=(const MatchboxPhasespace &) = delete;
361 
362 };
363 
364 }
365 
366 #endif /* HERWIG_MatchboxPhasespace_H */
const pair< size_t, size_t > & lastSingularIndices() const
Return the last matched singular limit.
bool useMassGenerators() const
True, if mass generators should be used instead of fixed masses.
Ptr< MatchboxPhasespace >::ptr cloneMe() const
Clone this phase space generator.
Energy singularCutoff
A cutoff below which a region is considered singular.
virtual bool haveX1X2() const
Return true, if this phasespace generator will generate incoming partons itself.
size_t nRnd
The number of random numbers available.
MatchboxPhasespace defines an abstract interface to a phase space generator.
virtual bool wantCMS() const
Return true, if this phase space generator expects the incoming partons in their center-of-mass syste...
double diagramWeight(const Tree2toNDiagram &diag) const
Return the weight appropriate to the given diagram.
virtual void setXComb(tStdXCombPtr xc)
Set the XComb object steering the Born matrix element this class represents virtual corrections to...
int theLoopParticleIdMin
The first id in a range of id&#39;s meant to denote fictitious &#39;ghost&#39; particles to be used by the diagra...
Provide easy access to MatchboxXComb XComb extensions.
void singularLimit(size_t i, size_t j)
Limit phasespace generation to a given collinear or soft limit.
Ptr< PhasespaceCouplings >::ptr theCouplings
Couplings to be used in diagram weighting.
void clearDiagramWeights()
Clear the diagram weights.
vector< cPDPtr > cPDVector
int theLoopParticleIdMax
The last id in a range of id&#39;s meant to denote fictitious &#39;ghost&#39; particles to be used by the diagram...
const double * numbers
The random numbers.
double operator()()
Return next random number.
vector< DiagPtr > DiagramVector
StreamingRnd()
Default constructor.
-*- C++ -*-
Wrap around a vector of random numbers to behave as a stream of those.
bool theUseMassGenerators
True, if mass generators should be used instead of fixed masses.
virtual double invertTwoToNKinematics(const vector< Lorentz5Momentum > &, double *) const
Invert the given phase space point to the random numbers which would have generated it...
virtual bool isInvertible() const
Return true, if this phase space generator is invertible.
StreamingRnd(const double *newNumbers, size_t n)
Construct from random numbers.