Page MenuHomeHEPForge
Files F19245073

jets.cc
No OneTemporary
Actions

Size
5 KB
Referenced Files
None
Subscribers
None

jets.cc
View Options

/**
* \authors The HEJ collaboration (see AUTHORS for details)
* \date 2019-2023
* \copyright GPLv2 or later
*/
#include "HEJ/jets.hh"
#include <algorithm>
#include <cassert>
#include <cmath>
#include "HEJ/Constants.hh"
// generated headers
#include "HEJ/currents/j_j.hh"
#include "HEJ/currents/jqqbar_j.hh"
#include "HEJ/currents/j_qqbar_j.hh"
#include "HEJ/currents/juno_j.hh"
namespace {
using COM = std::complex<double>;
} // Anonymous Namespace
namespace HEJ {
namespace currents {
double ME_qq(
HLV const & p1out, HLV const & p1in,
HLV const & p2out, HLV const & p2in
){
using helicity::plus;
using helicity::minus;
// we only have to distinguish between the two possibilities of
// contracting same-helicity or opposite-helicity currents.
COM const Mp = HEJ::j_j<minus, plus>(p1in, p1out, p2in, p2out);
COM const Mm = HEJ::j_j<minus, minus>(p1in, p1out, p2in, p2out);
double const sst = std::norm(Mm) + std::norm(Mp);
// Multiply by factor 2 (helicities)
return 2.*sst;
}
namespace {
template<Helicity h1, Helicity hg>
double amp_juno_j(
HLV const & pa, HLV const & pb,
HLV const & pg, HLV const & p1, HLV const & p2
) {
// TODO: code duplication with Wjets
const COM u1 = U1_j<h1, hg>(pa,p1,pb,p2,pg);
const COM u2 = U2_j<h1, hg>(pa,p1,pb,p2,pg);
const COM l = L_j<h1, hg>(pa,p1,pb,p2,pg);
const COM x = u1 - l;
const COM y = u2 + l;
return C_F*norm(x + y) - C_A*(x*conj(y)).real();
}
} // Anonymous Namespace
double ME_unob_qq(
HLV const & pg,
HLV const & p1out, HLV const & p1in,
HLV const & p2out, HLV const & p2in
){
using helicity::minus;
using helicity::plus;
// only calculate half of the helicity amplitudes,
// the remaining ones follow from CP symmetry
const double amm = amp_juno_j<minus, minus>(p1in, p2in, pg, p1out, p2out);
const double amp = amp_juno_j<minus, plus>(p1in, p2in, pg, p1out, p2out);
const double apm = amp_juno_j< plus, minus>(p1in, p2in, pg, p1out, p2out);
const double app = amp_juno_j< plus, plus>(p1in, p2in, pg, p1out, p2out);
constexpr double hel_fac = 2.;
return hel_fac * (amm + amp + apm + app);
}
namespace {
// helicity amplitudes for jqqbar j contraction
template<Helicity h1, Helicity hg>
double amp_jqqbar_j(
HLV const & pa, HLV const & pb,
HLV const & p1, HLV const & p2, HLV const & p3
) {
const COM u1 = U1X_j<h1, hg>(pa,pb,p1,p2,p3);
const COM u2 = U2X_j<h1, hg>(pa,pb,p1,p2,p3);
const COM l = LX_j <h1, hg>(pa,pb,p1,p2,p3);
const COM x = u1 - l;
const COM y = u2 + l;
return C_F*norm(x + y) - C_A*(x*conj(y)).real();
}
} // Anonymous Namespace
double ME_qqbar_qg(
HLV const & pa, HLV const & pb,
HLV const & p1, HLV const & p2, HLV const & p3
){
using helicity::minus;
using helicity::plus;
// only calculate half of the helicity amplitudes,
// the remaining ones follow from CP symmetry
const double ampsq =
+ amp_jqqbar_j<minus, minus>(pa, pb, p1, p2, p3)
+ amp_jqqbar_j<minus, plus> (pa, pb, p1, p2, p3)
+ amp_jqqbar_j<plus, minus>(pa, pb, p1, p2, p3)
+ amp_jqqbar_j<plus, plus> (pa, pb, p1, p2, p3)
;
constexpr double hel_fac = 2.;
// correct colour averaging after crossing:
constexpr double avg_fac = N_C / (N_C*N_C - 1.);
return avg_fac * hel_fac * ampsq;
}
namespace {
// helicity amplitudes for matrix element with central qqbar
template<Helicity h1a, Helicity hqqbar>
double amp_Cenqqbar_qq(
HLV const & pa, HLV const & p1,
HLV const & pb, HLV const & pn,
HLV const & pq, HLV const & pqbar,
HLV const & q11, HLV const & q12
){
using std::norm;
const COM Xs = M_Xs<h1a, hqqbar>(pa, p1, pb, pn, pq, pqbar, q11, q12);
const COM X6 = M_X6<h1a, hqqbar>(pa, p1, pb, pn, pq, pqbar, q11, q12);
const COM X7 = M_X7<h1a, hqqbar>(pa, p1, pb, pn, pq, pqbar, q11, q12);
const COM V = Xs + X6;
const COM W = Xs + X7;
return (C_A*C_F*C_F/4.) * (norm(V) + norm(W)) + (C_F/4.) * real(V * conj(W));
}
} // Anonymous Namespace
double ME_Cenqqbar_qq(
HLV const & pa, HLV const & pb,
std::vector<HLV> const & partons,
const bool qbar_first, const std::size_t nabove
){
using helicity::plus;
using helicity::minus;
HLV q1 = pa;
for(std::size_t i = 0; i <= nabove; ++i){
q1 -= partons[i];
}
const auto qq = split_into_lightlike(q1);
// since q1.m2() < 0 the following assertion is always true
// see documentation for split_into_lightlike
assert(qq.second.e() < 0);
const HLV pq = (qbar_first ? partons[nabove+2] : partons[nabove+1]);
const HLV pqbar = (qbar_first ? partons[nabove+1] : partons[nabove+2]);
const HLV p1 = partons.front();
const HLV pn = partons.back();
// 8 helicity choices, but only 4 independent ones
//(complex conjugation related).
// In principle, the proper helicity labeling depends on
// whether we have antiquark lines (aqlinea and aqlineb).
// However, this only corresponds to a relabeling,
// which doesn't change the sum over all helicities below.
const double ampsq =
+ amp_Cenqqbar_qq<minus, minus>(pa, p1, pb, pn, pq, pqbar, qq.first, -qq.second)
+ amp_Cenqqbar_qq<minus, plus> (pa, p1, pb, pn, pq, pqbar, qq.first, -qq.second)
+ amp_Cenqqbar_qq<plus, minus>(pa, p1, pb, pn, pq, pqbar, qq.first, -qq.second)
+ amp_Cenqqbar_qq<plus, plus> (pa, p1, pb, pn, pq, pqbar, qq.first, -qq.second);
//Result (averaged, without coupling or t-channel props). Factor of
//2 for the 4 helicity configurations I didn't work out explicitly
const HLV q3 = q1 - pq - pqbar;
constexpr double hel_fac = 2.;
// factor 1/2: see eq:Scen in developer manual
return hel_fac * ampsq / (2. * q1.m2() * q3.m2());
}
} // namespace currents
} // namespace HEJ

File Metadata

Mime Type
text/x-c
Expires
Tue, Sep 30, 4:48 AM (7 h, 11 m)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
6461238
Default Alt Text
jets.cc (5 KB)

Event Timeline

HEPForge.org · .