diff --git a/src/Wjets.cc b/src/Wjets.cc
index c63eabf..c72e6b8 100644
--- a/src/Wjets.cc
+++ b/src/Wjets.cc
@@ -1,1223 +1,1215 @@
/**
* \authors The HEJ collaboration (see AUTHORS for details)
* \date 2019
* \copyright GPLv2 or later
*/
#include "HEJ/jets.hh"
#include "HEJ/Wjets.hh"
#include "HEJ/Tensor.hh"
#include "HEJ/Constants.hh"
#include <array>
#include <iostream>
using HEJ::Tensor;
using HEJ::init_sigma_index;
using HEJ::metric;
using HEJ::rank3_current;
using HEJ::rank5_current;
using HEJ::eps;
using HEJ::to_tensor;
using HEJ::Helicity;
using HEJ::angle;
using HEJ::square;
namespace helicity = HEJ::helicity;
namespace { // Helper Functions
// FKL W Helper Functions
double WProp (const HLV & plbar, const HLV & pl){
COM propW = COM(0.,-1.)/( (pl+plbar).m2() - HEJ::MW*HEJ::MW + COM(0.,1.)*HEJ::MW*HEJ::GammaW);
double PropFactor=(propW*conj(propW)).real();
return PropFactor;
}
CCurrent jW (
HLV pout, Helicity helout,
HLV plbar, Helicity hellbar,
HLV pl, Helicity hell,
HLV pin, Helicity helin
){
COM cur[4];
cur[0]=0.;
cur[1]=0.;
cur[2]=0.;
cur[3]=0.;
CCurrent sum(0.,0.,0.,0.);
// NOTA BENE: Conventions for W+ --> e+ nu, so that nu is lepton(6), e is
// anti-lepton(5)
// Need to swap e and nu for events with W- --> e- nubar!
if (helin==helout && hellbar==hell) {
HLV qa=pout+plbar+pl;
HLV qb=pin-plbar-pl;
double ta(qa.m2()),tb(qb.m2());
CCurrent temp2,temp3,temp5;
CCurrent t65 = joo(pl,hell,plbar,hellbar);
CCurrent vout(pout.e(),pout.x(),pout.y(),pout.z());
CCurrent vin(pin.e(),pin.x(),pin.y(),pin.z());
COM brac615=t65.dot(vout);
COM brac645=t65.dot(vin);
// prod1565 and prod6465 are zero for Ws (not Zs)!!
temp2 = joo(pout,helout,pl,helout);
COM prod1665=temp2.dot(t65);
temp3 = joi(plbar,helin,pin,helin);
COM prod5465=temp3.dot(t65);
temp2=joo(pout,helout,plbar,helout);
temp3=joi(pl,hell,pin,helin);
temp5=joi(pout,helout,pin,helin);
CCurrent term1,term2,term3;
term1=(2.*brac615/ta+2.*brac645/tb)*temp5;
term2=(prod1665/ta)*temp3;
term3=(-prod5465/tb)*temp2;
sum=term1+term2+term3;
}
return sum;
}
namespace {
// see eq. (87) {eq:jW-} in developer manual
Tensor<1> jW_minus(
HLV const & pa, HLV const & p1,
HLV const & plbar, HLV const & pl
){
using HEJ::helicity::minus;
const double tWin = (pa-pl-plbar).m2();
const double tWout = (p1+pl+plbar).m2();
constexpr COM i{0,1};
// C++ arithmetic operators are evaluated left-to-right,
// so the following first computes complex scalar coefficients,
// which then multiply a current, reducing the number
// of multiplications
return 2.*i*(
+ angle(p1, pl)*square(p1, plbar)/tWout
+ square(pa, plbar)*angle(pa, pl)/tWin
)*HEJ::current(p1, pa, helicity::minus)
+ 2.*i*angle(p1, pl)*square(pl, plbar)/tWout
*HEJ::current(pl, pa, helicity::minus)
+ 2.*i*square(pa, plbar)*angle(pl, plbar)/tWin
*HEJ::current(p1, plbar, helicity::minus);
}
}
// see eqs. (87), (88) {eq:jW-}, {eq:jW+} in developer manual
Tensor<1> jW(
HLV const & pa, HLV const & p1,
HLV const & plbar, HLV const & pl,
Helicity h
){
if(h == helicity::minus) {
return jW_minus(pa, p1, plbar, pl);
}
return jW_minus(pa, p1, pl, plbar).complex_conj();
}
CCurrent jWbar (
HLV pout, Helicity helout,
HLV plbar, Helicity hellbar,
HLV pl, Helicity hell,
HLV pin, Helicity helin
){
COM cur[4];
cur[0]=0.;
cur[1]=0.;
cur[2]=0.;
cur[3]=0.;
CCurrent sum(0.,0.,0.,0.);
// NOTA BENE: Conventions for W+ --> e+ nu, so that nu is lepton(6), e is
// anti-lepton(5)
// Need to swap e and nu for events with W- --> e- nubar!
if (helin==helout && hellbar==hell) {
HLV qa=pout+plbar+pl;
HLV qb=pin-plbar-pl;
double ta(qa.m2()),tb(qb.m2());
CCurrent temp2,temp3,temp5;
CCurrent t65 = joo(pl,hell,plbar,hellbar);
CCurrent vout(pout.e(),pout.x(),pout.y(),pout.z());
CCurrent vin(pin.e(),pin.x(),pin.y(),pin.z());
COM brac615=t65.dot(vout);
COM brac645=t65.dot(vin);
// prod1565 and prod6465 are zero for Ws (not Zs)!!
temp2 = joo(plbar,helout,pout,helout); // temp2 is <5|alpha|1>
COM prod5165=temp2.dot(t65);
temp3 = jio(pin,helin,pl,helin); // temp3 is <4|alpha|6>
COM prod4665=temp3.dot(t65);
temp2=joo(pl,helout,pout,helout); // temp2 is now <6|mu|1>
temp3=jio(pin,helin,plbar,helin); // temp3 is now <4|mu|5>
temp5=jio(pin,helin,pout,helout); // temp5 is <4|mu|1>
CCurrent term1,term2,term3;
term1 =(-2.*brac615/ta-2.*brac645/tb)*temp5;
term2 =(-prod5165/ta)*temp3;
term3 =(prod4665/tb)*temp2;
sum = term1 + term2 + term3;
}
return sum;
}
/**
* @brief W+Jets Unordered Contribution Helper Functions
* @returns result of equation (4.1.28) in Helen's Thesis (p.100)
*/
double jM2Wuno(HLV pg, HLV p1,HLV plbar, HLV pl, HLV pa, Helicity h1,
HLV p2, HLV pb, Helicity h2, Helicity pol
){
//@TODO Simplify the below (less Tensor class?)
init_sigma_index();
HLV pW = pl+plbar;
HLV q1g=pa-pW-p1-pg;
HLV q1 = pa-p1-pW;
HLV q2 = p2-pb;
const double taW = (pa-pW).m2();
const double taW1 = (pa-pW-p1).m2();
const double tb2 = (pb-p2).m2();
const double tb2g = (pb-p2-pg).m2();
const double s1W = (p1+pW).m2();
const double s1gW = (p1+pW+pg).m2();
const double s1g = (p1+pg).m2();
const double tag = (pa-pg).m2();
const double taWg = (pa-pW-pg).m2();
//use p1 as ref vec in pol tensor
Tensor<1> epsg = eps(pg,p2,pol);
Tensor<1> epsW = HEJ::current(pl,plbar,helicity::minus);
Tensor<1> j2b = HEJ::current(p2,pb,h2);
Tensor<1> Tq1q2 = to_tensor((q1+q2)/taW1 + (pb/pb.dot(pg)
+p2/p2.dot(pg)) * tb2/(2*tb2g));
Tensor<3> J31a = rank3_current(p1, pa, h1);
Tensor<2> J2_qaW =J31a.contract((pa-pW)/taW, 2);
Tensor<2> J2_p1W =J31a.contract((p1+pW)/s1W, 2);
Tensor<3> L1a = outer(Tq1q2, J2_qaW);
Tensor<3> L1b = outer(Tq1q2, J2_p1W);
Tensor<3> L2a = outer(-pg-q1,J2_qaW)/taW1;
Tensor<3> L2b = outer(-pg-q1, J2_p1W)/taW1;
Tensor<3> L3 = outer(metric(), J2_qaW.contract(pg-q2,1)+J2_p1W.contract(pg-q2,2))/taW1;
Tensor<3> L(0.);
Tensor<5> J51a = rank5_current(p1, pa, h1);
Tensor<4> J_qaW = J51a.contract((pa-pW),4);
Tensor<4> J_qag = J51a.contract(pa-pg,4);
Tensor<4> J_p1gW = J51a.contract(p1+pg+pW,4);
Tensor<3> U1a = J_qaW.contract(p1+pg,2);
Tensor<3> U1b = J_p1gW.contract(p1+pg,2);
Tensor<3> U1c = J_p1gW.contract(p1+pW,2);
Tensor<3> U1(0.);
Tensor<3> U2a = J_qaW.contract(pa-pg-pW,2);
Tensor<3> U2b = J_qag.contract(pa-pg-pW,2);
Tensor<3> U2c = J_qag.contract(p1+pW,2);
Tensor<3> U2(0.);
for(int nu=0; nu<4;nu++){
for(int mu=0;mu<4;mu++){
for(int rho=0;rho<4;rho++){
L(nu, mu, rho) = L1a(nu,mu,rho) + L1b(nu,rho,mu)
+ L2a(mu,nu,rho) + L2b(mu,rho,nu) + L3(mu,nu,rho);
U1(nu, mu, rho) = U1a(nu, mu, rho) / (s1g*taW)
+ U1b(nu,rho,mu) / (s1g*s1gW) + U1c(rho,nu,mu) / (s1W*s1gW);
U2(nu,mu,rho) = U2a(mu,nu,rho) / (taWg*taW)
+ U2b(mu,rho,nu) / (taWg*tag) + U2c(rho,mu,nu) / (s1W*tag);
}
}
}
COM X = ((((U1-L).contract(epsW,3)).contract(j2b,2)).contract(epsg,1));
COM Y = ((((U2+L).contract(epsW,3)).contract(j2b,2)).contract(epsg,1));
double amp = HEJ::C_A*HEJ::C_F*HEJ::C_F/2.*(norm(X)+norm(Y)) - HEJ::C_F/2.*(X*conj(Y)).real();
double t1 = q1g.m2();
double t2 = q2.m2();
double WPropfact = WProp(plbar, pl);
//Divide by WProp
amp*=WPropfact;
//Divide by t-channels
amp/=(t1*t2);
return amp;
}
// Relevant Wqqx Helper Functions.
//g->qxqlxl (Calculates gluon to qqx Current. See JV_\mu in WSubleading Notes)
Tensor <1> gtqqxW(HLV pq,HLV pqbar,HLV pl,HLV plbar){
//@TODO Simplify the calculation below (Less Tensor class use?)
double s2AB=(pl+plbar+pq).m2();
double s3AB=(pl+plbar+pqbar).m2();
// Define llx current.
Tensor<1> ABCur = HEJ::current(pl, plbar, helicity::minus);
//blank 3 Gamma Current
Tensor<3> JV23 = rank3_current(pq,pqbar,helicity::minus);
// Components of g->qqW before W Contraction
Tensor<2> JV1 = JV23.contract((pq + pl + plbar),2)/(s2AB);
Tensor<2> JV2 = JV23.contract((pqbar + pl + plbar),2)/(s3AB);
// g->qqW Current. Note Minus between terms due to momentum flow.
// Also note: (-I)^2 from W vert. (I) from Quark prop.
Tensor<1> JVCur = (JV1.contract(ABCur,1) - JV2.contract(ABCur,2))*COM(0.,-1.);
return JVCur;
}
// Helper Functions Calculate the Crossed Contribution
Tensor <2> MCrossW(HLV pa, HLV, HLV, HLV, HLV pq, HLV pqbar, HLV pl,
HLV plbar, std::vector<HLV> partons, int nabove
){
//@TODO Simplify the calculation below Maybe combine with MCross?
// Useful propagator factors
double s2AB=(pl+plbar+pq).m2();
double s3AB=(pl+plbar+pqbar).m2();
HLV q1, q3;
q1=pa;
for(int i=0; i<nabove+1;i++){
q1=q1-partons.at(i);
}
q3 = q1 - pq - pqbar - pl - plbar;
double tcro1=(q3+pq).m2();
double tcro2=(q1-pqbar).m2();
// Define llx current.
Tensor<1> ABCur = HEJ::current(pl, plbar,helicity::minus);
//Blank 5 gamma Current
Tensor<5> J523 = rank5_current(pq,pqbar,helicity::minus);
// 4 gamma currents (with 1 contraction already).
Tensor<4> J_q3q = J523.contract((q3 + pq),2);
Tensor<4> J_2AB = J523.contract((pq + pl + plbar),2);
// Components of Crossed Vertex Contribution
Tensor<3> Xcro1 = J_q3q.contract((pqbar + pl + plbar),3);
Tensor<3> Xcro2 = J_q3q.contract((q1 - pqbar),3);
Tensor<3> Xcro3 = J_2AB.contract((q1 - pqbar),3);
// Term Denominators Taken Care of at this stage
Tensor<2> Xcro1Cont = Xcro1.contract(ABCur,3)/(tcro1*s3AB);
Tensor<2> Xcro2Cont = Xcro2.contract(ABCur,2)/(tcro1*tcro2);
Tensor<2> Xcro3Cont = Xcro3.contract(ABCur,1)/(s2AB*tcro2);
//Initialise the Crossed Vertex Object
Tensor<2> Xcro(0.);
for(int mu=0; mu<4;mu++){
for(int nu=0;nu<4;nu++){
Xcro(mu,nu) = -(-Xcro1Cont(nu,mu)+Xcro2Cont(nu,mu)+Xcro3Cont(nu,mu));
}
}
return Xcro;
}
// Helper Functions Calculate the Uncrossed Contribution
Tensor <2> MUncrossW(HLV pa, HLV, HLV, HLV, HLV pq, HLV pqbar,
HLV pl, HLV plbar, std::vector<HLV> partons, int nabove
){
//@TODO Simplify the calculation below Maybe combine with MUncross?
double s2AB=(pl+plbar+pq).m2();
double s3AB=(pl+plbar+pqbar).m2();
HLV q1, q3;
q1=pa;
for(int i=0; i<nabove+1;i++){
q1=q1-partons.at(i);
}
q3 = q1 - pl - plbar - pq - pqbar;
double tunc1 = (q1-pq).m2();
double tunc2 = (q3+pqbar).m2();
// Define llx current.
Tensor<1> ABCur = HEJ::current(pl, plbar, helicity::minus);
//Blank 5 gamma Current
Tensor<5> J523 = rank5_current(pq,pqbar,helicity::minus);
// 4 gamma currents (with 1 contraction already).
Tensor<4> J_2AB = J523.contract((pq + pl + plbar),2);
Tensor<4> J_q1q = J523.contract((q1 - pq),2);
// 2 Contractions taken care of.
Tensor<3> Xunc1 = J_2AB.contract((q3 + pqbar),3);
Tensor<3> Xunc2 = J_q1q.contract((q3 + pqbar),3);
Tensor<3> Xunc3 = J_q1q.contract((pqbar + pl + plbar),3);
// Term Denominators Taken Care of at this stage
Tensor<2> Xunc1Cont = Xunc1.contract(ABCur,1)/(s2AB*tunc2);
Tensor<2> Xunc2Cont = Xunc2.contract(ABCur,2)/(tunc1*tunc2);
Tensor<2> Xunc3Cont = Xunc3.contract(ABCur,3)/(tunc1*s3AB);
//Initialise the Uncrossed Vertex Object
Tensor<2> Xunc(0.);
for(int mu=0; mu<4;mu++){
for(int nu=0;nu<4;nu++){
Xunc(mu,nu) = -(- Xunc1Cont(mu,nu)+Xunc2Cont(mu,nu) +Xunc3Cont(mu,nu));
}
}
return Xunc;
}
// Helper Functions Calculate the g->qqxW (Eikonal) Contributions
Tensor <2> MSymW(HLV pa, HLV p1, HLV pb, HLV p4, HLV pq, HLV pqbar,
HLV pl,HLV plbar, std::vector<HLV> partons, int nabove
){
//@TODO Simplify the calculation below Maybe combine with MSym?
double sa2=(pa+pq).m2();
double s12=(p1+pq).m2();
double sa3=(pa+pqbar).m2();
double s13=(p1+pqbar).m2();
double saA=(pa+pl).m2();
double s1A=(p1+pl).m2();
double saB=(pa+plbar).m2();
double s1B=(p1+plbar).m2();
double sb2=(pb+pq).m2();
double s42=(p4+pq).m2();
double sb3=(pb+pqbar).m2();
double s43=(p4+pqbar).m2();
double sbA=(pb+pl).m2();
double s4A=(p4+pl).m2();
double sbB=(pb+plbar).m2();
double s4B=(p4+plbar).m2();
double s23AB=(pl+plbar+pq+pqbar).m2();
HLV q1,q3;
q1=pa;
for(int i=0;i<nabove+1;i++){
q1-=partons.at(i);
}
q3=q1-pq-pqbar-plbar-pl;
double t1 = (q1).m2();
double t3 = (q3).m2();
// g->qqW Current (Factors of sqrt2 dealt with in this function.)
Tensor<1> JV = gtqqxW(pq,pqbar,pl,plbar);
// 1a gluon emisson Contribution
Tensor<3> X1a = outer(metric(), p1*(t1/(s12+s13+s1A+s1B))
+ pa*(t1/(sa2+sa3+saA+saB)) );
Tensor<2> X1aCont = X1a.contract(JV,3);
//4b gluon emission Contribution
Tensor<3> X4b = outer(metric(), p4*(t3/(s42+s43+s4A+s4B))
+ pb*(t3/(sb2+sb3+sbA+sbB)) );
Tensor<2> X4bCont = X4b.contract(JV,3);
//Set up each term of 3G diagram.
Tensor<3> X3g1 = outer(q1+pq+pqbar+pl+plbar, metric());
Tensor<3> X3g2 = outer(q3-pq-pqbar-pl-plbar, metric());
Tensor<3> X3g3 = outer(q1+q3, metric());
// Note the contraction of indices changes term by term
Tensor<2> X3g1Cont = X3g1.contract(JV,3);
Tensor<2> X3g2Cont = X3g2.contract(JV,2);
Tensor<2> X3g3Cont = X3g3.contract(JV,1);
// XSym is an amalgamation of x1a, X4b and X3g.
// Makes sense from a colour factor point of view.
Tensor<2>Xsym(0.);
for(int mu=0; mu<4;mu++){
for(int nu=0;nu<4;nu++){
Xsym(mu,nu) = (X3g1Cont(nu,mu) + X3g2Cont(mu,nu) - X3g3Cont(nu,mu))
+ (X1aCont(mu,nu) - X4bCont(mu,nu));
}
}
return Xsym/s23AB;
}
Tensor <2> MCross(HLV pa, HLV pq, HLV pqbar, std::vector<HLV> partons,
Helicity hq, int nabove
){
//@TODO Simplify the calculation below Maybe combine with MCrossW?
HLV q1;
q1=pa;
for(int i=0;i<nabove+1;i++){
q1-=partons.at(i);
}
double t2=(q1-pqbar).m2();
//Blank 3 gamma Current
Tensor<3> J323 = rank3_current(pq,pqbar,hq);
// 2 gamma current (with 1 contraction already).
Tensor<2> XCroCont = J323.contract((q1-pqbar),2)/(t2);
//Initialise the Crossed Vertex
Tensor<2> Xcro(0.);
for(int mu=0; mu<4;mu++){
for(int nu=0;nu<4;nu++){
Xcro(mu,nu) = XCroCont(nu,mu);
}
}
return Xcro;
}
// Helper Functions Calculate the Uncrossed Contribution
Tensor <2> MUncross(HLV pa, HLV pq,HLV pqbar, std::vector<HLV> partons,
Helicity hq, int nabove
){
//@TODO Simplify the calculation below Maybe combine with MUncrossW?
HLV q1;
q1=pa;
for(int i=0;i<nabove+1;i++){
q1-=partons.at(i);
}
double t2 = (q1-pq).m2();
//Blank 3 gamma Current
Tensor<3> J323 = rank3_current(pq,pqbar,hq);
// 2 gamma currents (with 1 contraction already).
Tensor<2> XUncCont = J323.contract((q1-pq),2)/t2;
//Initialise the Uncrossed Vertex
Tensor<2> Xunc(0.);
for(int mu=0; mu<4;mu++){
for(int nu=0;nu<4;nu++){
Xunc(mu,nu) = -XUncCont(mu,nu);
}
}
return Xunc;
}
// Helper Functions Calculate the Eikonal Contributions
Tensor <2> MSym(HLV pa, HLV p1, HLV pb, HLV p4, HLV pq, HLV pqbar,
std::vector<HLV> partons, Helicity hq, int nabove
){
//@TODO Simplify the calculation below Maybe combine with MsymW?
HLV q1, q3;
q1=pa;
for(int i=0;i<nabove+1;i++){
q1-=partons.at(i);
}
q3 = q1-pq-pqbar;
double t1 = (q1).m2();
double t3 = (q3).m2();
double s23 = (pq+pqbar).m2();
double sa2 = (pa+pq).m2();
double sa3 = (pa+pqbar).m2();
double s12 = (p1+pq).m2();
double s13 = (p1+pqbar).m2();
double sb2 = (pb+pq).m2();
double sb3 = (pb+pqbar).m2();
double s42 = (p4+pq).m2();
double s43 = (p4+pqbar).m2();
Tensor<1> qqxCur = HEJ::current(pq, pqbar, hq);
// // 1a gluon emisson Contribution
Tensor<3> X1a = outer(metric(), p1*(t1/(s12+s13))+ pa*(t1/(sa2+sa3)));
Tensor<2> X1aCont = X1a.contract(qqxCur,3);
// //4b gluon emission Contribution
Tensor<3> X4b = outer(metric(), p4*(t3/(s42+s43)) + pb*(t3/(sb2+sb3)));
Tensor<2> X4bCont = X4b.contract(qqxCur,3);
// New Formulation Corresponding to New Analytics
Tensor<3> X3g1 = outer(q1+pq+pqbar, metric());
Tensor<3> X3g2 = outer(q3-pq-pqbar, metric());
Tensor<3> X3g3 = outer(q1+q3, metric());
// Note the contraction of indices changes term by term
Tensor<2> X3g1Cont = X3g1.contract(qqxCur,3);
Tensor<2> X3g2Cont = X3g2.contract(qqxCur,2);
Tensor<2> X3g3Cont = X3g3.contract(qqxCur,1);
Tensor<2>Xsym(0.);
for(int mu=0; mu<4;mu++){
for(int nu=0;nu<4;nu++){
Xsym(mu, nu) = COM(0,1) * ( (X3g1Cont(nu,mu) + X3g2Cont(mu,nu)
- X3g3Cont(nu,mu)) + (X1aCont(mu,nu) - X4bCont(mu,nu)) );
}
}
return Xsym/s23;
}
//! W+Jets FKL Contributions
/**
* @brief W+Jets FKL Contributions, function to handle all incoming types.
* @param p1out Outgoing Particle 1. (W emission)
* @param plbar Outgoing election momenta
* @param pl Outgoing neutrino momenta
* @param p1in Incoming Particle 1. (W emission)
* @param p2out Outgoing Particle 2
* @param p2in Incoming Particle 2
* @param aqlineb Bool. Is Backwards quark line an anti-quark line?
* @param aqlinef Bool. Is Forwards quark line an anti-quark line?
*
* Calculates j_W ^\mu j_\mu.
* Handles all possible incoming states.
*/
double jW_j( HLV p1out, HLV plbar, HLV pl, HLV p1in, HLV p2out, HLV p2in,
- bool aqlineb, bool aqlinef
+ bool aqlineb, bool /* aqlinef */
){
- CCurrent mj1m,mj2p,mj2m;
- HLV q1=p1in-p1out-plbar-pl;
- HLV q2=-(p2in-p2out);
-
- if(aqlineb) mj1m=jWbar(p1out,helicity::minus, plbar,helicity::minus,
- pl,helicity::minus, p1in,helicity::minus);
- else mj1m=jW(p1out,helicity::minus, plbar,helicity::minus,
- pl,helicity::minus, p1in,helicity::minus);
-
- mj2p=joi(p2out,!aqlinef, p2in,!aqlinef);
- mj2m=joi(p2out, aqlinef, p2in, aqlinef);
-
- COM Mmp=mj1m.dot(mj2p);
- COM Mmm=mj1m.dot(mj2m);
-
- // sum of spinor strings ||^2
- double a2Mmp=abs2(Mmp);
- double a2Mmm=abs2(Mmm);
-
- double WPropfact = WProp(plbar, pl);
-
+ using helicity::minus;
+ using helicity::plus;
+ const HLV q1=p1in-p1out-plbar-pl;
+ const HLV q2=-(p2in-p2out);
+
+ const double WPropfact = WProp(plbar, pl);
+
+ const auto j_W = COM{0,-1}*jW(p1in, p1out, plbar, pl, aqlineb?plus:minus);
+ double Msqr = 0.;
+ for(const auto h: {plus, minus}) {
+ const auto j = HEJ::current(p2out, p2in, h);
+ Msqr += abs2(j_W.contract(j, 1));
+ }
// Division by colour and Helicity average (Nc2-1)(4)
// Multiply by Cf^2
- return HEJ::C_F*HEJ::C_F*WPropfact*(a2Mmp+a2Mmm)/(q1.m2()*q2.m2()*(HEJ::N_C*HEJ::N_C - 1)*4);
+ return HEJ::C_F*HEJ::C_F*WPropfact*Msqr/(q1.m2()*q2.m2()*(HEJ::N_C*HEJ::N_C - 1)*4);
}
} // Anonymous Namespace
double ME_W_qQ (HLV p1out, HLV plbar, HLV pl,HLV p1in, HLV p2out, HLV p2in){
return jW_j(p1out, plbar, pl, p1in, p2out, p2in, false, false);
}
double ME_W_qQbar (HLV p1out, HLV plbar, HLV pl,HLV p1in, HLV p2out, HLV p2in){
return jW_j(p1out, plbar, pl, p1in, p2out, p2in, false, true);
}
double ME_W_qbarQ (HLV p1out, HLV plbar, HLV pl,HLV p1in, HLV p2out, HLV p2in){
return jW_j(p1out, plbar, pl, p1in, p2out, p2in, true, false);
}
double ME_W_qbarQbar (HLV p1out, HLV plbar, HLV pl,HLV p1in, HLV p2out, HLV p2in){
return jW_j(p1out, plbar, pl, p1in, p2out, p2in, true, true);
}
double ME_W_qg (HLV p1out, HLV plbar, HLV pl,HLV p1in, HLV p2out, HLV p2in){
return jW_j(p1out, plbar, pl, p1in, p2out, p2in, false, false)*K_g(p2out, p2in)/HEJ::C_F;
}
double ME_W_qbarg (HLV p1out, HLV plbar, HLV pl,HLV p1in, HLV p2out, HLV p2in){
return jW_j(p1out, plbar, pl, p1in, p2out, p2in, true, false)*K_g(p2out, p2in)/HEJ::C_F;
}
namespace{
/**
* @brief W+Jets Unordered Contributions, function to handle all incoming types.
* @param p1out Outgoing Particle 1. (W emission)
* @param plbar Outgoing election momenta
* @param pl Outgoing neutrino momenta
* @param p1in Incoming Particle 1. (W emission)
* @param p2out Outgoing Particle 2 (Quark, unordered emission this side.)
* @param p2in Incoming Particle 2 (Quark, unordered emission this side.)
* @param pg Unordered Gluon momenta
* @param aqlineb Bool. Is Backwards quark line an anti-quark line?
* @param aqlinef Bool. Is Forwards quark line an anti-quark line?
*
* Calculates j_W ^\mu j_{uno}_\mu. Ie, unordered with W emission opposite side.
* Handles all possible incoming states.
*/
double jW_juno(HLV p1out, HLV plbar, HLV pl,HLV p1in, HLV p2out,
HLV p2in, HLV pg, bool aqlineb, bool aqlinef){
CCurrent mj1m,mj2p,mj2m, jgbm,jgbp,j2gm,j2gp;
HLV q1=p1in-p1out-plbar-pl;
HLV q2=-(p2in-p2out-pg);
HLV q3=-(p2in-p2out);
if(aqlineb) mj1m=jWbar(p1out,helicity::minus, plbar,helicity::minus,
pl,helicity::minus, p1in,helicity::minus);
else mj1m=jW(p1out,helicity::minus, plbar,helicity::minus,
pl,helicity::minus, p1in,helicity::minus);
// Note <p1|eps|pa> current split into two by gauge choice.
// See James C's Thesis (p72). <p1|eps|pa> -> <p1|pg><pg|pa>
mj2p=joi(p2out,!aqlinef,p2in,!aqlinef);
mj2m=joi(p2out,aqlinef,p2in,aqlinef);
jgbp=joi(pg,!aqlinef,p2in,!aqlinef);
jgbm=joi(pg,aqlinef,p2in,aqlinef);
// Note for function joo(): <p1+|pg+> = <pg-|p1->.
j2gp=joo(p2out,!aqlinef,pg,!aqlinef);
j2gm=joo(p2out,aqlinef,pg,aqlinef);
// Dot products of these which occur again and again
COM MWmp=mj1m.dot(mj2p); // And now for the Higgs ones
COM MWmm=mj1m.dot(mj2m);
CCurrent qsum(q2+q3);
CCurrent Lmp,Lmm,Lpp,Lpm,U1mp,U1mm,U1pp,U1pm,U2mp,U2mm,U2pp,U2pm,p1o(p1out),p1i(p1in);
CCurrent p2o(p2out);
CCurrent p2i(p2in);
Lmm=( (-1.)*qsum*(MWmm) + (-2.*mj1m.dot(pg))*mj2m + 2.*mj2m.dot(pg)*mj1m
+ ( p1o/pg.dot(p1out) + p1i/pg.dot(p1in) )*( q2.m2()*MWmm/2. ) )/q3.m2();
Lmp=( (-1.)*qsum*(MWmp) + (-2.*mj1m.dot(pg))*mj2p + 2.*mj2p.dot(pg)*mj1m
+ ( p1o/pg.dot(p1out) + p1i/pg.dot(p1in) )*( q2.m2()*MWmp/2. ) )/q3.m2();
U1mm=(jgbm.dot(mj1m)*j2gm+2.*p2o*MWmm)/(p2out+pg).m2();
U1mp=(jgbp.dot(mj1m)*j2gp+2.*p2o*MWmp)/(p2out+pg).m2();
U2mm=((-1.)*j2gm.dot(mj1m)*jgbm+2.*p2i*MWmm)/(p2in-pg).m2();
U2mp=((-1.)*j2gp.dot(mj1m)*jgbp+2.*p2i*MWmp)/(p2in-pg).m2();
double amm,amp;
amm=HEJ::C_F*(2.*vre(Lmm-U1mm,Lmm+U2mm))+2.*HEJ::C_F*HEJ::C_F/3.*vabs2(U1mm+U2mm);
amp=HEJ::C_F*(2.*vre(Lmp-U1mp,Lmp+U2mp))+2.*HEJ::C_F*HEJ::C_F/3.*vabs2(U1mp+U2mp);
double ampsq=-(amm+amp);
//Divide by WProp
ampsq*=WProp(plbar, pl);
return ampsq/((16)*(q2.m2()*q1.m2()));
}
}
double ME_W_unob_qQ(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jW_juno(p2out, plbar, pl, p2in, p1out, p1in, pg, false, false);
}
double ME_W_unob_qQbar(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jW_juno(p2out, plbar, pl, p2in, p1out, p1in, pg, false, true);
}
double ME_W_unob_qbarQ(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jW_juno(p2out, plbar, pl, p2in, p1out, p1in, pg, true, false);
}
double ME_W_unob_qbarQbar(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jW_juno(p2out, plbar, pl, p2in, p1out, p1in, pg, true, true);
}
double ME_W_unof_qQ(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jW_juno(p1out, plbar, pl, p1in, p2out, p2in, pg, false, false);
}
double ME_W_unof_qQbar(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jW_juno(p1out, plbar, pl, p1in, p2out, p2in, pg, false, true);
}
double ME_W_unof_qbarQ(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jW_juno(p1out, plbar, pl, p1in, p2out, p2in, pg, true, false);
}
double ME_W_unof_qbarQbar(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jW_juno(p1out, plbar, pl, p1in, p2out, p2in, pg, true, true);
}
namespace{
/**
* @brief W+Jets Unordered Contributions, function to handle all incoming types.
* @param pg Unordered Gluon momenta
* @param p1out Outgoing Particle 1. (Quark - W and Uno emission)
* @param plbar Outgoing election momenta
* @param pl Outgoing neutrino momenta
* @param p1in Incoming Particle 1. (Quark - W and Uno emission)
* @param p2out Outgoing Particle 2
* @param p2in Incoming Particle 2
* @param aqlineb Bool. Is Backwards quark line an anti-quark line?
*
* Calculates j_W_{uno} ^\mu j_\mu. Ie, unordered with W emission same side.
* Handles all possible incoming states. Note this handles both forward and back-
* -ward Wuno emission. For forward, ensure p1out is the uno and W emission parton.
* @TODO: Include separate wrapper functions for forward and backward to clean up
* ME_W_unof_current in `MatrixElement.cc`.
*/
double jWuno_j(HLV pg, HLV p1out, HLV plbar, HLV pl, HLV p1in,
HLV p2out, HLV p2in, bool aqlineb
){
//Calculate different Helicity choices
const Helicity h = aqlineb?helicity::plus:helicity::minus;
double ME2mpp = jM2Wuno(pg, p1out,plbar,pl,p1in,h,p2out,p2in,helicity::plus,helicity::plus);
double ME2mpm = jM2Wuno(pg, p1out,plbar,pl,p1in,h,p2out,p2in,helicity::plus,helicity::minus);
double ME2mmp = jM2Wuno(pg, p1out,plbar,pl,p1in,h,p2out,p2in,helicity::minus,helicity::plus);
double ME2mmm = jM2Wuno(pg, p1out,plbar,pl,p1in,h,p2out,p2in,helicity::minus,helicity::minus);
//Helicity sum and average over initial states
return (ME2mpp + ME2mpm + ME2mmp + ME2mmm)/(4.*HEJ::C_A*HEJ::C_A);
}
}
double ME_Wuno_qQ(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jWuno_j(pg, p1out, plbar, pl, p1in, p2out, p2in, false);
}
double ME_Wuno_qQbar(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jWuno_j(pg, p1out, plbar, pl, p1in, p2out, p2in, false);
}
double ME_Wuno_qbarQ(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jWuno_j(pg, p1out, plbar, pl, p1in, p2out, p2in, true);
}
double ME_Wuno_qbarQbar(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jWuno_j(pg, p1out, plbar, pl, p1in, p2out, p2in, true);
}
double ME_Wuno_qg(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jWuno_j(pg, p1out, plbar, pl, p1in, p2out, p2in, false)*K_g(p2out, p2in)/HEJ::C_F;
}
double ME_Wuno_qbarg(HLV p1out, HLV p1in, HLV p2out, HLV p2in,
HLV pg, HLV plbar, HLV pl
){
return jWuno_j(pg, p1out, plbar, pl, p1in, p2out, p2in, true)*K_g(p2out, p2in)/HEJ::C_F;
}
/**
* @brief W+Jets Extremal qqx Contributions, function to handle all incoming types.
* @param pgin Incoming gluon which will split into qqx.
* @param pqout Quark of extremal qqx outgoing (W-Emission).
* @param plbar Outgoing anti-lepton momenta
* @param pl Outgoing lepton momenta
* @param pqbarout Anti-quark of extremal qqx pair. (W-Emission)
* @param pout Outgoing Particle 2 (end of FKL chain)
* @param p2in Incoming Particle 2
* @param aqlinef Bool. Is Forwards quark line an anti-quark line?
*
* Calculates j_W_{qqx} ^\mu j_\mu. Ie, Ex-QQX with W emission same side.
* Handles all possible incoming states. Calculated via crossing symmetry from jWuno_j.
*/
double jWqqx_j(HLV pgin, HLV pqout, HLV plbar, HLV pl,
HLV pqbarout, HLV p2out, HLV p2in, bool aqlinef){
//Calculate Different Helicity Configurations.
const Helicity h = aqlinef?helicity::plus:helicity::minus;
double ME2mpp = jM2Wuno(-pgin, pqout,plbar,pl,-pqbarout,h,p2out,p2in,helicity::plus,helicity::plus);
double ME2mpm = jM2Wuno(-pgin, pqout,plbar,pl,-pqbarout,h,p2out,p2in,helicity::plus,helicity::minus);
double ME2mmp = jM2Wuno(-pgin, pqout,plbar,pl,-pqbarout,h,p2out,p2in,helicity::minus,helicity::plus);
double ME2mmm = jM2Wuno(-pgin, pqout,plbar,pl,-pqbarout,h,p2out,p2in,helicity::minus,helicity::minus);
//Helicity sum and average over initial states.
double ME2 = (ME2mpp + ME2mpm + ME2mmp + ME2mmm)/(4.*HEJ::C_A*HEJ::C_A);
//Correct colour averaging after crossing:
ME2*=(3.0/8.0);
return ME2;
}
double ME_WExqqx_qbarqQ(HLV pgin, HLV pqout, HLV plbar, HLV pl,
HLV pqbarout, HLV p2out, HLV p2in){
return jWqqx_j(pgin, pqout, plbar, pl, pqbarout, p2out, p2in, false);
}
double ME_WExqqx_qqbarQ(HLV pgin, HLV pqbarout, HLV plbar, HLV pl,
HLV pqout, HLV p2out, HLV p2in){
return jWqqx_j(pgin, pqbarout, plbar, pl, pqout, p2out, p2in, true);
}
double ME_WExqqx_qbarqg(HLV pgin, HLV pqout, HLV plbar, HLV pl,
HLV pqbarout, HLV p2out, HLV p2in){
return jWqqx_j(pgin, pqout, plbar, pl, pqbarout, p2out, p2in, false)*K_g(p2out,p2in)/HEJ::C_F;
}
double ME_WExqqx_qqbarg(HLV pgin, HLV pqbarout, HLV plbar, HLV pl,
HLV pqout, HLV p2out, HLV p2in){
return jWqqx_j(pgin, pqbarout, plbar, pl, pqout, p2out, p2in, true)*K_g(p2out,p2in)/HEJ::C_F;
}
namespace {
//Function to calculate Term 1 in Equation 3.23 in James Cockburn's Thesis.
Tensor<1> qggm1(HLV pb, HLV p2, HLV p3, Helicity hel2, Helicity helg, HLV refmom){
//@TODO Simplify the calculation below. (Less Tensor class use?)
double t1 = (p3-pb)*(p3-pb);
// Gauge choice in polarisation tensor. (see JC's Thesis)
Tensor<1> epsg = eps(pb, refmom, helg);
Tensor<3> qqCurBlank = rank3_current(p2,p3,hel2);
Tensor<2> qqCur = qqCurBlank.contract(p3-pb,2);
Tensor<1> gqqCur = qqCur.contract(epsg,2)/t1;
return gqqCur*(-1);
}
//Function to calculate Term 2 in Equation 3.23 in James Cockburn's Thesis.
Tensor<1> qggm2(HLV pb, HLV p2, HLV p3, Helicity hel2, Helicity helg, HLV refmom){
//@TODO Simplify the calculation below (Less Tensor class use?)
double t1 = (p2-pb)*(p2-pb);
// Gauge choice in polarisation tensor. (see JC's Thesis)
Tensor<1> epsg = eps(pb,refmom, helg);
Tensor<3> qqCurBlank = rank3_current(p2,p3,hel2);
Tensor<2> qqCur = qqCurBlank.contract(p2-pb,2);
Tensor<1> gqqCur = qqCur.contract(epsg,1)/t1;
return gqqCur;
}
//Function to calculate Term 3 in Equation 3.23 in James Cockburn's Thesis.
Tensor<1> qggm3(HLV pb, HLV p2, HLV p3, Helicity hel2, Helicity helg, HLV refmom){
//@TODO Simplify the calculation below (Less Tensor class use?)
double s23 = (p2+p3)*(p2+p3);
// Gauge choice in polarisation tensor. (see JC's Thesis)
Tensor<1> epsg = eps(pb, refmom, helg);
Tensor<3> qqCurBlank1 = outer(p2+p3, metric())/s23;
Tensor<3> qqCurBlank2 = outer(pb, metric())/s23;
Tensor<1> Cur23 = HEJ::current(p2, p3,hel2);
Tensor<2> qqCur1 = qqCurBlank1.contract(Cur23,3);
Tensor<2> qqCur2 = qqCurBlank2.contract(Cur23,3);
Tensor<2> qqCur3 = qqCurBlank2.contract(Cur23,1);
Tensor<1> gqqCur = (qqCur1.contract(epsg,1)
- qqCur2.contract(epsg,2)
+ qqCur3.contract(epsg,1))*2*COM(0,1);
return gqqCur;
}
}
// no wqq emission
double ME_W_Exqqx_QQq(HLV pa, HLV pb, HLV p1, HLV p2,
HLV p3,HLV plbar, HLV pl, bool aqlinepa
){
using helicity::minus;
using helicity::plus;
init_sigma_index();
// 2 independent helicity choices (complex conjugation related).
Tensor<1> TMmmm1 = qggm1(pb,p2,p3,minus,minus, pa);
Tensor<1> TMmmm2 = qggm2(pb,p2,p3,minus,minus, pa);
Tensor<1> TMmmm3 = qggm3(pb,p2,p3,minus,minus, pa);
Tensor<1> TMpmm1 = qggm1(pb,p2,p3,minus,plus, pa);
Tensor<1> TMpmm2 = qggm2(pb,p2,p3,minus,plus, pa);
Tensor<1> TMpmm3 = qggm3(pb,p2,p3,minus,plus, pa);
// Build the external quark line W Emmision
Tensor<1> cur1a = jW(pa,p1,plbar,pl, aqlinepa?plus:minus);
//Contract with the qqxCurrent.
COM Mmmm1 = TMmmm1.contract(cur1a,1);
COM Mmmm2 = TMmmm2.contract(cur1a,1);
COM Mmmm3 = TMmmm3.contract(cur1a,1);
COM Mpmm1 = TMpmm1.contract(cur1a,1);
COM Mpmm2 = TMpmm2.contract(cur1a,1);
COM Mpmm3 = TMpmm3.contract(cur1a,1);
//Colour factors:
COM cm1m1,cm2m2,cm3m3,cm1m2,cm1m3,cm2m3;
cm1m1=8./3.;
cm2m2=8./3.;
cm3m3=6.;
cm1m2 =-1./3.;
cm1m3 = -3.*COM(0.,1.);
cm2m3 = 3.*COM(0.,1.);
//Sqaure and sum for each helicity config:
double Mmmm = real( cm1m1*pow(abs(Mmmm1),2) + cm2m2*pow(abs(Mmmm2),2)
+ cm3m3*pow(abs(Mmmm3),2) + 2.*real(cm1m2*Mmmm1*conj(Mmmm2))
+ 2.*real(cm1m3*Mmmm1*conj(Mmmm3))
+ 2.*real(cm2m3*Mmmm2*conj(Mmmm3)) );
double Mpmm = real( cm1m1*pow(abs(Mpmm1),2) + cm2m2*pow(abs(Mpmm2),2)
+ cm3m3*pow(abs(Mpmm3),2) + 2.*real(cm1m2*Mpmm1*conj(Mpmm2))
+ 2.*real(cm1m3*Mpmm1*conj(Mpmm3))
+ 2.*real(cm2m3*Mpmm2*conj(Mpmm3)) );
// Divide by WProp
double WPropfact = WProp(plbar, pl);
return (2*WPropfact*(Mmmm+Mpmm)/24./4.)/(pa-p1-pl-plbar).m2()/(p2+p3-pb).m2();
}
// W+Jets qqxCentral
double ME_WCenqqx_qq(HLV pa, HLV pb,HLV pl, HLV plbar, std::vector<HLV> partons,
bool aqlinepa, bool aqlinepb, bool qqxmarker, int nabove
){
init_sigma_index();
HLV pq, pqbar, p1, p4;
if (qqxmarker){
pqbar = partons[nabove+1];
pq = partons[nabove+2];}
else{
pq = partons[nabove+1];
pqbar = partons[nabove+2];}
p1 = partons.front();
p4 = partons.back();
Tensor<1> T1am, T4bm, T1ap, T4bp;
if(!(aqlinepa)){
T1ap = HEJ::current(p1, pa, helicity::plus);
T1am = HEJ::current(p1, pa, helicity::minus);}
else if(aqlinepa){
T1ap = HEJ::current(pa, p1, helicity::plus);
T1am = HEJ::current(pa, p1, helicity::minus);}
if(!(aqlinepb)){
T4bp = HEJ::current(p4, pb, helicity::plus);
T4bm = HEJ::current(p4, pb, helicity::minus);}
else if(aqlinepb){
T4bp = HEJ::current(pb, p4, helicity::plus);
T4bm = HEJ::current(pb, p4, helicity::minus);}
// Calculate the 3 separate contributions to the effective vertex
Tensor<2> Xunc = MUncrossW(pa, p1, pb, p4, pq, pqbar, pl, plbar, partons, nabove);
Tensor<2> Xcro = MCrossW( pa, p1, pb, p4, pq, pqbar, pl, plbar, partons, nabove);
Tensor<2> Xsym = MSymW( pa, p1, pb, p4, pq, pqbar, pl, plbar, partons, nabove);
// 4 Different Helicity Choices (Differs from Pure Jet Case, where there is
// also the choice in qqbar helicity.
// (- - hel choice)
COM M_mmUnc = (((Xunc).contract(T1am,1)).contract(T4bm,1));
COM M_mmCro = (((Xcro).contract(T1am,1)).contract(T4bm,1));
COM M_mmSym = (((Xsym).contract(T1am,1)).contract(T4bm,1));
// (- + hel choice)
COM M_mpUnc = (((Xunc).contract(T1am,1)).contract(T4bp,1));
COM M_mpCro = (((Xcro).contract(T1am,1)).contract(T4bp,1));
COM M_mpSym = (((Xsym).contract(T1am,1)).contract(T4bp,1));
// (+ - hel choice)
COM M_pmUnc = (((Xunc).contract(T1ap,1)).contract(T4bm,1));
COM M_pmCro = (((Xcro).contract(T1ap,1)).contract(T4bm,1));
COM M_pmSym = (((Xsym).contract(T1ap,1)).contract(T4bm,1));
// (+ + hel choice)
COM M_ppUnc = (((Xunc).contract(T1ap,1)).contract(T4bp,1));
COM M_ppCro = (((Xcro).contract(T1ap,1)).contract(T4bp,1));
COM M_ppSym = (((Xsym).contract(T1ap,1)).contract(T4bp,1));
//Colour factors:
COM cmsms,cmumu,cmcmc,cmsmu,cmsmc,cmumc;
cmsms=3.;
cmumu=4./3.;
cmcmc=4./3.;
cmsmu =3./2.*COM(0.,1.);
cmsmc = -3./2.*COM(0.,1.);
cmumc = -1./6.;
// Work Out Interference in each case of helicity:
double amp_mm = real(cmsms*pow(abs(M_mmSym),2)
+cmumu*pow(abs(M_mmUnc),2)
+cmcmc*pow(abs(M_mmCro),2)
+2.*real(cmsmu*M_mmSym*conj(M_mmUnc))
+2.*real(cmsmc*M_mmSym*conj(M_mmCro))
+2.*real(cmumc*M_mmUnc*conj(M_mmCro)));
double amp_mp = real(cmsms*pow(abs(M_mpSym),2)
+cmumu*pow(abs(M_mpUnc),2)
+cmcmc*pow(abs(M_mpCro),2)
+2.*real(cmsmu*M_mpSym*conj(M_mpUnc))
+2.*real(cmsmc*M_mpSym*conj(M_mpCro))
+2.*real(cmumc*M_mpUnc*conj(M_mpCro)));
double amp_pm = real(cmsms*pow(abs(M_pmSym),2)
+cmumu*pow(abs(M_pmUnc),2)
+cmcmc*pow(abs(M_pmCro),2)
+2.*real(cmsmu*M_pmSym*conj(M_pmUnc))
+2.*real(cmsmc*M_pmSym*conj(M_pmCro))
+2.*real(cmumc*M_pmUnc*conj(M_pmCro)));
double amp_pp = real(cmsms*pow(abs(M_ppSym),2)
+cmumu*pow(abs(M_ppUnc),2)
+cmcmc*pow(abs(M_ppCro),2)
+2.*real(cmsmu*M_ppSym*conj(M_ppUnc))
+2.*real(cmsmc*M_ppSym*conj(M_ppCro))
+2.*real(cmumc*M_ppUnc*conj(M_ppCro)));
double amp=((amp_mm+amp_mp+amp_pm+amp_pp)/(9.*4.));
HLV q1,q3;
q1=pa;
for(int i=0;i<nabove+1;i++){
q1-=partons.at(i);
}
q3 = q1 - pq - pqbar - pl - plbar;
double t1 = (q1).m2();
double t3 = (q3).m2();
//Divide by t-channels
amp/=(t1*t1*t3*t3);
//Divide by WProp
double WPropfact = WProp(plbar, pl);
amp*=WPropfact;
return amp;
}
// no wqq emission
double ME_W_Cenqqx_qq(HLV pa, HLV pb,HLV pl,HLV plbar, std::vector<HLV> partons,
bool aqlinepa, bool aqlinepb, bool qqxmarker, int nabove,
int nbelow, bool forwards
){
using helicity::minus;
using helicity::plus;
init_sigma_index();
if (!forwards){ //If Emission from Leg a instead, flip process.
std::swap(pa, pb);
std::reverse(partons.begin(),partons.end());
std::swap(aqlinepa, aqlinepb);
qqxmarker = !qqxmarker;
std::swap(nabove,nbelow);
}
HLV pq, pqbar, p1,p4;
if (qqxmarker){
pqbar = partons[nabove+1];
pq = partons[nabove+2];}
else{
pq = partons[nabove+1];
pqbar = partons[nabove+2];}
p1 = partons.front();
p4 = partons.back();
Tensor<1> T1am(0.), T1ap(0.);
if(!(aqlinepa)){
T1ap = HEJ::current(p1, pa, plus);
T1am = HEJ::current(p1, pa, minus);}
else if(aqlinepa){
T1ap = HEJ::current(pa, p1, plus);
T1am = HEJ::current(pa, p1, minus);}
Tensor <1> T4bm = jW(pb, p4, plbar, pl, aqlinepb?plus:minus);
// Calculate the 3 separate contributions to the effective vertex
Tensor<2> Xunc_m = MUncross(pa, pq, pqbar,partons, minus, nabove);
Tensor<2> Xcro_m = MCross( pa, pq, pqbar,partons, minus, nabove);
Tensor<2> Xsym_m = MSym( pa, p1, pb, p4, pq, pqbar, partons, minus, nabove);
Tensor<2> Xunc_p = MUncross(pa, pq, pqbar,partons, plus, nabove);
Tensor<2> Xcro_p = MCross( pa, pq, pqbar,partons, plus, nabove);
Tensor<2> Xsym_p = MSym( pa, p1, pb, p4, pq, pqbar, partons, plus, nabove);
// (- - hel choice)
COM M_mmUnc = (((Xunc_m).contract(T1am,1)).contract(T4bm,1));
COM M_mmCro = (((Xcro_m).contract(T1am,1)).contract(T4bm,1));
COM M_mmSym = (((Xsym_m).contract(T1am,1)).contract(T4bm,1));
// (- + hel choice)
COM M_mpUnc = (((Xunc_p).contract(T1am,1)).contract(T4bm,1));
COM M_mpCro = (((Xcro_p).contract(T1am,1)).contract(T4bm,1));
COM M_mpSym = (((Xsym_p).contract(T1am,1)).contract(T4bm,1));
// (+ - hel choice)
COM M_pmUnc = (((Xunc_m).contract(T1ap,1)).contract(T4bm,1));
COM M_pmCro = (((Xcro_m).contract(T1ap,1)).contract(T4bm,1));
COM M_pmSym = (((Xsym_m).contract(T1ap,1)).contract(T4bm,1));
// (+ + hel choice)
COM M_ppUnc = (((Xunc_p).contract(T1ap,1)).contract(T4bm,1));
COM M_ppCro = (((Xcro_p).contract(T1ap,1)).contract(T4bm,1));
COM M_ppSym = (((Xsym_p).contract(T1ap,1)).contract(T4bm,1));
//Colour factors:
COM cmsms,cmumu,cmcmc,cmsmu,cmsmc,cmumc;
cmsms=3.;
cmumu=4./3.;
cmcmc=4./3.;
cmsmu =3./2.*COM(0.,1.);
cmsmc = -3./2.*COM(0.,1.);
cmumc = -1./6.;
// Work Out Interference in each case of helicity:
double amp_mm = real(cmsms*pow(abs(M_mmSym),2)
+cmumu*pow(abs(M_mmUnc),2)
+cmcmc*pow(abs(M_mmCro),2)
+2.*real(cmsmu*M_mmSym*conj(M_mmUnc))
+2.*real(cmsmc*M_mmSym*conj(M_mmCro))
+2.*real(cmumc*M_mmUnc*conj(M_mmCro)));
double amp_mp = real(cmsms*pow(abs(M_mpSym),2)
+cmumu*pow(abs(M_mpUnc),2)
+cmcmc*pow(abs(M_mpCro),2)
+2.*real(cmsmu*M_mpSym*conj(M_mpUnc))
+2.*real(cmsmc*M_mpSym*conj(M_mpCro))
+2.*real(cmumc*M_mpUnc*conj(M_mpCro)));
double amp_pm = real(cmsms*pow(abs(M_pmSym),2)
+cmumu*pow(abs(M_pmUnc),2)
+cmcmc*pow(abs(M_pmCro),2)
+2.*real(cmsmu*M_pmSym*conj(M_pmUnc))
+2.*real(cmsmc*M_pmSym*conj(M_pmCro))
+2.*real(cmumc*M_pmUnc*conj(M_pmCro)));
double amp_pp = real(cmsms*pow(abs(M_ppSym),2)
+cmumu*pow(abs(M_ppUnc),2)
+cmcmc*pow(abs(M_ppCro),2)
+2.*real(cmsmu*M_ppSym*conj(M_ppUnc))
+2.*real(cmsmc*M_ppSym*conj(M_ppCro))
+2.*real(cmumc*M_ppUnc*conj(M_ppCro)));
double amp=((amp_mm+amp_mp+amp_pm+amp_pp)/(9.*4.));
HLV q1,q3;
q1=pa;
for(int i=0;i<nabove+1;i++){
q1-=partons.at(i);
}
q3 = q1 - pq - pqbar;
double t1 = (q1).m2();
double t3 = (q3).m2();
//Divide by t-channels
amp/=(t1*t1*t3*t3);
//Divide by WProp
double WPropfact = WProp(plbar, pl);
amp*=WPropfact;
return amp;
}
diff --git a/t/ME_data/ME_Wm.dat b/t/ME_data/ME_Wm.dat
index fffffa4..5c554d1 100644
--- a/t/ME_data/ME_Wm.dat
+++ b/t/ME_data/ME_Wm.dat
@@ -1,1308 +1,1308 @@
2.928675692e-06
1.647116461e-06
6.792853649e-07
8.309213874e-05
1.188586835e-09
7.418850985e-09
4.641858351e-07
4.469299635e-07
6.109075964e-06
2.937428029e-06
4.834759305e-06
8.223080997e-06
1.984882889e-06
4.724373663e-06
2.284512955e-07
6.086791509e-08
9.402879575e-07
2.134932663e-05
1.344153838e-06
6.325190057e-08
5.042205822e-06
0.0002285013435
2.294857015e-06
6.277759204e-05
1.337020968e-07
3.953835541e-06
9.658506379e-05
2.316868337e-06
2.819688016e-07
0.0001081216619
3.432998402e-06
3.368523741e-06
1.11288709e-08
3.29543466e-07
3.989570061e-10
1.948796723e-05
3.510830208e-07
3.210202477e-07
3.239728865e-06
8.604875852e-06
7.617842534e-10
1.878886612e-06
8.602552313e-09
2.905134545e-09
9.580831652e-08
2.761377426e-07
1.121758088e-05
8.305715929e-07
3.094575366e-06
9.009051399e-07
3.660116555e-07
4.545389317e-08
1.619299497e-09
0.01367636241
3.915774786e-06
2.025205952e-06
3.6814996e-09
4.179895596e-06
3.07283947e-08
8.89627591e-07
3.025959937e-08
1.817403251e-09
4.337439968e-10
1.271896091e-06
6.045011555e-06
9.89332115e-09
1.214106087e-05
2.958683161e-08
2.061281825e-09
1.93033348e-05
6.299313702e-07
2.644016282e-07
2.061367847e-06
4.34253878e-05
8.291732575e-07
2.499007952e-06
4.358489891e-06
9.241699824e-05
9.12437674e-05
0.0001159481718
2.350252878e-06
7.845324212e-06
0.000638536879
4.723440426e-08
6.954486577e-09
5.357700283e-06
1.473730225e-07
1.054004829e-07
5.920067571e-07
1.175811475e-05
1.546221661e-07
0.0004236100189
3.617420289e-06
4.516879429e-05
3.955718614e-05
0.0001365328248
5.126170478e-05
1.077742001e-06
1.018758481e-08
2.310396111e-09
6.836738619e-07
5.827603078e-09
1.557658774e-09
6.048954972e-05
1.652858793e-07
2.956993046e-07
3.715591231e-06
1.403931806e-08
2.060713606e-07
0.001033146046
9.773118638e-07
2.935985027e-05
1.151256054e-05
7.794713809e-06
2.17436307e-05
4.056751186e-08
2.955377882e-08
8.030872004e-05
2.603215792e-07
0.0001681672442
2.930133644e-08
7.197036394e-05
1.920989813e-06
3.637607439e-09
7.893842383e-06
0.0001175217074
4.802925073e-10
4.83556567e-07
2.878827354e-05
1.898723086e-09
2.754516812e-08
1.403679092e-10
5.826558306e-07
2.433773671e-06
7.927708105e-06
1.020983626e-08
1.052585256e-07
7.737192672e-06
3.400645867e-06
2.714553382e-08
2.377796617e-09
0.0001490868071
9.612832588e-11
1.654742249e-05
1.57073488e-07
1.609407803e-06
6.946367688e-08
0.04337790658
9.860795665e-06
1.103405546e-05
0.01558805851
1.359140493e-06
2.104988229e-07
4.483328013e-08
2.527261731e-08
4.136266286e-07
1.75962763e-07
4.16236019e-07
3.66072855e-07
2.475097776e-06
2.142069091e-05
3.962234552e-07
5.0752339e-06
2.831972358e-08
8.260228591e-08
2.191466894e-08
7.392657036e-05
5.57330681e-10
1.264513056e-08
3.373173459e-05
1.515657759e-07
4.183645033e-07
1.162296215e-07
7.731815775e-06
1.077601875e-05
2.331663061e-07
4.316286429e-10
4.05516744e-09
1.282346148e-09
1.905147472e-09
2.799718364e-07
6.430044079e-08
3.336162544e-07
2.847516226e-07
2.8800296e-05
5.034816791e-08
0.0001697859258
7.357868068e-08
1.312972542e-07
1.494668933e-07
1.309682993e-09
1.452620265e-08
7.003475182e-07
2.392867531e-06
6.228302547e-06
2.32207831e-07
8.172158542e-08
4.289261389e-07
2.619009165e-07
7.684975193e-06
3.9300806e-10
1.235480059e-06
4.345702611e-07
1.019591398e-07
1.014400081e-05
1.198511446e-05
3.072393212e-05
2.47145635e-06
6.141459748e-06
2.217791841e-05
7.506816715e-07
4.156597855e-07
0.0003928655193
6.192512495e-07
9.058620171e-07
7.767884384e-06
1.506541078e-06
1.438366886e-08
1.393919618e-07
7.338405876e-05
1.013840918e-07
1.941618603e-05
1.098271271e-06
7.583485063e-07
1.924040495e-08
2.683414048e-06
0.0005387206665
4.514003014e-09
7.345041607e-06
4.863384343e-08
5.137710227e-07
1.341025387e-07
1.385096075e-09
0.001475731456
4.766899706e-09
2.788152153e-07
4.885059085e-07
4.389966638e-07
0.0006056993831
2.260959244e-06
1.622050932e-07
8.044376305e-07
0.002874337098
8.545204208e-05
2.859915875e-07
3.520051428e-05
3.83840617e-05
6.701394914e-12
6.711555934e-09
1.221714284e-07
4.175432874e-09
1.276099272e-06
9.349373741e-06
5.093349519e-07
2.102059997e-05
2.546943624e-06
4.307372366e-07
0.0001405081607
1.532593492e-06
0.004327081818
1.798381994e-07
0.0002232914245
1.177488117e-07
1.784201554e-06
2.027515818e-09
3.74919864e-07
1.128537898e-07
4.135413677e-05
2.009104635e-08
8.773550962e-07
9.824802351e-10
1.631528339e-07
2.288354101e-05
2.035248598e-07
0.0001019103655
3.912289372e-06
5.329003459e-08
3.815264848e-10
1.388142926e-05
2.066904312e-08
9.173930152e-07
6.957230348e-09
2.903877008e-05
4.478794954e-09
6.770693006e-07
9.714079428e-08
1.998883986e-05
0.001799109057
1.236217788e-08
1.118375624e-05
1.144499288e-07
1.1864756e-06
4.716258378e-06
3.464456734e-07
7.532453996e-11
1.742294659e-06
1.345426643e-06
8.449545648e-08
8.954221651e-07
2.432487759e-06
0.0001400677698
0.001713444386
1.075048353e-07
4.105406541e-05
3.733656311e-07
1.984952182e-06
0.001480322688
0.0001151217974
0.00209287627
3.27555809e-06
2.185782336e-09
0.0004548139237
0.01477027961
3.934273672e-05
0.001243829892
5.430764382e-08
4.380856265e-07
4.297458002e-06
1.467632945e-08
1.820453663e-08
2.699514835e-09
8.420814307e-07
8.634913084e-10
1.182745993e-06
8.508146264e-09
7.184712053e-07
1.192590185e-06
1.766408443e-08
0.002934313619
1.283384372e-06
6.255055172e-06
8.703331507e-06
9.664774229e-06
2.675288802e-09
3.392650797e-08
0.0001071682868
5.832833131e-09
8.881117172e-08
0.0002166070842
1.248608916e-08
1.726759584e-06
2.99705455e-06
2.119854874e-07
7.412646571e-06
2.755998704e-08
9.960376424e-06
3.051416349e-05
1.129492506e-08
2.619064731e-08
5.532504512e-10
2.177158595e-10
0.0003851956976
4.247896087e-05
1.005180535e-08
7.063067217e-08
0.0002453874122
6.569452237e-09
2.430119913e-06
1.727937098e-07
1.211730869e-05
9.94795514e-06
0.0001605757277
5.961378248e-08
1.63764285e-06
1.925056588e-07
1.895751671e-05
3.584139499e-09
2.424148932e-07
0.001006931825
1.003162717e-06
1.237484464e-06
2.006432918e-06
1.665677334e-05
8.731387573e-06
1.333871332e-07
7.301399559e-09
0.0005093528313
7.645422963e-08
2.020202554e-06
6.813091876e-08
1.10814416e-06
1.95551597e-05
2.003879505e-06
3.670933614e-07
1.728622505e-07
0.0007263797041
-1.411925493e-05
+1.41191122028095e-05
1.571311147e-07
3.392287113e-06
1.377047714e-09
6.630650805e-07
5.227904356e-08
1.127011492e-06
0.0003754771874
1.865169444e-07
6.576952888e-06
3.329497044e-06
7.07257004e-06
2.88154104e-07
1.985441412e-07
6.929088895e-08
2.512260458e-06
6.217791051e-06
3.877741985e-07
1.707903155e-08
5.10233078e-08
4.572653519e-09
0.1095009268
7.449132821e-07
2.09740422e-06
9.378911822e-09
7.910046337e-08
1.223499404e-05
2.672193111e-10
2.314556691e-08
3.503343355e-13
1.602894349e-07
5.890632996e-07
1.487885964e-05
1.978704988e-05
2.552121299e-08
3.565602472e-08
0.0004467518296
9.443687582e-10
1.363026913e-06
2.235788016e-07
1.136525751e-07
2.962478329e-08
9.303967807e-06
3.925082515e-10
1.017880604e-05
1.434347997e-08
4.47759051e-06
1.996980966e-09
4.66863509e-05
5.845204451e-06
1.509519246e-07
1.741949098e-06
2.215116471e-10
1.532616641e-06
1.571651059e-08
7.107477006e-05
7.234687759e-06
1.621643207e-06
2.56453036e-08
8.069391642e-09
2.035933164e-07
-7.550817486e-08
+7.550717369310519e-08
8.659572881e-06
-0.000432153513
+0.0004321439671563004
3.337724504e-05
4.782082889e-05
2.323101105e-08
1.19309915e-05
0.0001350492201
2.75397447e-06
1.823126697e-08
2.411767675e-05
2.283494193e-09
3.766390067e-08
3.588343942e-05
6.787383948e-05
8.832490704e-10
7.778718617e-10
1.917361443e-08
5.849351718e-06
2.28541183e-06
6.146818441e-06
1.827359426e-05
3.231234221e-05
9.921186251e-07
3.779404373e-09
1.924149804e-05
7.698409214e-06
1.221205448e-06
1.642479046e-06
3.301215266e-07
0.0004602920172
3.293196018e-07
3.950685019e-07
2.805273751e-05
1.739207752e-08
3.773037585e-06
1.537008391e-06
2.434913209e-08
8.315956868e-05
6.88859927e-06
6.869368159e-08
4.672444205e-09
7.02639721e-05
7.305366536e-07
2.162269825e-05
0.0004328624246
1.653264992e-07
1.553411281e-05
5.165593047e-10
2.153731042e-10
8.592408471e-09
8.974979512e-09
7.520620911e-07
2.080333069e-07
1.653991504e-07
9.9043617e-09
6.542874961e-09
0.0005318553158
1.144555776e-05
6.936210932e-05
7.62010085e-08
1.834443898e-08
8.099206646e-06
1.946663573e-09
4.916072586e-06
2.013222915e-05
8.510967278e-07
3.442296018e-07
1.194853419e-07
1.251233575e-07
2.847893723e-06
2.519238034e-06
3.662564708e-07
1.646271735e-07
5.373793619e-07
4.353476735e-08
8.727968448e-08
1.053026857e-05
1.083603462e-06
1.461412043e-06
1.220147207e-06
0.0009307060856
3.48570189e-07
2.208034886e-08
8.740398335e-08
0.0001005182181
5.658935843e-06
3.660460092e-07
5.210218655e-06
1.176510494e-09
1.167510222e-09
1.74723538e-07
4.652418777e-10
1.639422267e-06
9.200100093e-06
0.0001764752281
7.379070775e-09
2.542944745e-07
-3.71783565e-06
+3.717769041983323e-06
4.778249403e-10
6.09916388e-07
7.389397856e-06
1.030208854e-08
5.091853226e-08
1.929321992e-06
2.953703561e-07
1.428266933e-08
3.65646162e-05
8.358213962e-06
1.269052515e-05
3.142034606e-10
1.192354439e-08
3.42519028e-07
1.238073279e-05
4.592188983e-07
4.907760576e-08
2.605803561e-07
1.290266109e-06
8.91015308e-08
1.171147313e-05
6.948329894e-06
0.0001320627046
4.467361346e-08
1.343011103e-06
4.69742523e-09
1.674877668e-07
1.158457021e-07
7.247035325e-08
0.0001328992567
1.456012145e-05
1.974800428e-05
9.916997919e-05
7.076016655e-07
7.36330363e-07
8.750739052e-07
2.098702267e-05
1.224885671e-06
6.974041336e-07
6.433390899e-08
1.654868124e-05
0.0009302881479
1.821579155e-07
6.800514105e-07
1.766155755e-07
8.004650647e-07
9.941901599e-05
2.06752057e-05
5.066538235e-11
1.506772094e-06
8.10845767e-05
7.610845142e-05
0.0002115710954
2.12571763e-09
2.910738954e-06
9.48875277e-06
1.847390948e-06
7.436554496e-05
1.286439322e-10
1.413106449e-06
4.724231895e-08
5.001854028e-07
1.14055774e-06
6.718621271e-11
8.472867506e-09
8.494930984e-08
2.498427647e-07
6.987204698e-06
4.08117914e-06
2.82229796e-07
3.98589612e-09
6.278231145e-05
9.793860499e-07
9.666736058e-07
2.708628413e-06
3.490044251e-07
3.782645161e-10
8.131830282e-08
6.664325458e-06
2.365503287e-08
1.096958037e-06
3.776083109e-07
2.990539686e-08
0.002754346975
4.047644861e-07
1.411988389e-07
3.694067085e-08
0.0002320708655
3.71338597e-07
3.835998288e-05
2.623922776e-05
6.474281304e-10
0.0003219022096
3.356011806e-07
5.603569089e-09
9.856721257e-08
2.793861685e-07
5.195375958e-06
2.632243695e-06
6.350302361e-07
3.302377879e-10
3.043526423e-08
2.372803041e-08
6.617140209e-06
6.330476168e-07
1.870130608e-06
2.338555597e-07
4.591171378e-10
5.61933523e-05
1.955958771e-08
3.044175709e-08
5.161106748e-08
0.0001378680549
2.750711746e-09
2.289389014e-08
4.001113699e-09
2.236415233e-10
3.473425904e-06
7.311516826e-06
8.898096465e-08
1.481216022e-08
2.435208033e-12
1.189787819e-10
5.333724185e-06
4.664375968e-11
2.381386364e-10
1.249165957e-08
1.128630934e-07
1.550941521e-06
2.297705942e-12
5.846728434e-10
4.069354108e-10
4.063144762e-08
1.749554102e-08
5.155167146e-11
1.688297139e-11
2.774930205e-09
6.091667068e-08
2.87170322e-10
3.66643664e-10
8.258926512e-10
7.390910771e-09
1.48715158e-10
6.565786352e-12
3.475089984e-09
2.002380042e-09
7.274321744e-11
6.359749104e-09
7.940393308e-13
1.098611945e-06
7.438494044e-07
2.339207348e-10
4.499165346e-09
3.082882191e-09
1.878464872e-10
4.005585989e-09
2.63828057e-10
2.218789497e-09
5.827641219e-09
1.543332934e-10
4.907581247e-10
6.901620609e-10
2.088068542e-13
7.837109261e-10
1.937213353e-10
1.715044299e-08
4.069918057e-06
2.91951045e-08
1.565114165e-08
3.443423804e-08
1.256948998e-07
4.363190898e-09
2.287393528e-12
6.966944807e-09
5.995890076e-08
1.289835029e-10
1.695089305e-10
2.301282525e-11
9.445474271e-09
2.115405965e-08
1.662358545e-09
1.025068565e-09
9.719373728e-11
1.732767627e-11
4.41946888e-11
5.657095809e-11
-5.288646625e-07
+5.288528777290672e-07
5.324238867e-07
2.221589565e-07
2.48971658e-09
2.994534368049623e-09
1.075626423e-09
7.782686578e-11
2.946128307e-11
3.399399058e-11
1.446405065e-08
2.187225844e-08
7.820238975e-09
2.40229777e-09
7.46594076e-10
1.281010307e-06
3.598219943e-08
1.251224651967984e-08
3.332304338e-11
7.617797708e-10
9.909367591e-11
4.653906139e-10
3.587191007e-06
4.446718734e-08
1.894297187e-09
4.744577059e-05
4.433959517e-09
2.210894416e-08
1.036168503e-09
1.103345617e-07
4.995240386e-08
8.363740584e-12
3.281022886e-10
3.148002148e-08
6.330438798e-10
3.863775656e-10
5.96035067e-10
4.255627067e-11
1.165981634e-07
2.50543766e-09
1.33188621e-10
9.172189534e-07
6.582275731e-08
1.542558412e-10
1.297572348e-09
8.107196167e-09
2.070977659e-06
6.93677041e-09
2.329665548e-07
2.684341075e-07
3.375521064e-12
4.781691571e-11
5.350983865e-09
9.177711402e-10
3.059268527e-08
1.113405482e-07
5.998651665e-10
3.911876462e-10
1.176693635e-08
8.693402036e-10
2.320788866e-05
2.036112191e-07
3.21298594e-08
7.129850446e-08
1.481220144e-09
9.13398584e-07
9.416502413e-11
1.481517933e-08
2.317808347e-06
3.901257454e-08
4.627625989e-12
2.494572686e-08
9.094749026e-10
1.731947437e-08
2.095822035e-08
7.215954213e-08
2.372708431e-09
9.374174236e-13
1.029823717e-11
-5.170443463e-07
+5.170363884761981e-07
3.317414045e-07
1.109713675e-08
6.47170925e-07
1.864432618e-06
3.397272389e-11
4.168205664e-08
1.346830208e-11
3.171866066e-06
7.866352681e-11
3.289251567e-10
7.057976269e-14
2.67944287e-07
4.480095364e-10
2.542532019e-10
1.046877898e-09
4.692932633e-12
5.841384928e-10
5.254681082e-07
9.003591534e-11
3.637362188e-13
9.047013259e-09
1.769883631e-08
9.48770427e-07
1.037626188e-10
7.931261407e-08
2.511617984e-07
1.749603857e-09
1.158718175e-07
1.833652906e-06
4.450176914e-08
5.71537314e-10
1.070556403e-10
2.448737924e-10
6.262272951e-06
2.141707476e-08
4.258981326e-09
5.49246536e-09
9.566266675e-12
3.797543387e-11
2.705369876e-08
4.238910221e-07
1.784747449e-09
3.595414385e-08
3.222996111e-09
6.210164067e-05
2.258039015e-08
1.191514021e-13
1.513303142e-06
9.658853509e-11
5.082215114e-06
9.44022394e-10
8.437053744e-09
1.463048848e-10
3.444930056e-09
7.868562747e-12
1.272388225e-11
3.427574116e-10
8.990110201e-08
2.83369365e-11
3.0757007e-10
5.265605819e-14
4.26291208e-13
7.393316606e-13
5.385884913e-05
8.367399832e-10
1.271720599e-05
8.40960929e-11
6.150536572e-08
5.764038353e-10
1.4827533e-05
2.322813811e-10
4.431637476e-09
1.319799794e-08
1.965682244e-09
5.782112632e-06
7.526132679e-07
1.075618822e-09
1.712525448e-09
3.952010247e-07
3.136663182e-10
1.366784094e-10
2.291961377e-13
7.468862634e-11
2.13931093e-07
2.325391866e-09
7.657293509e-07
3.023206769e-07
1.974192114e-11
2.239591754e-13
1.575764712e-08
-3.894663093e-07
+3.894706527845038e-07
2.463427596e-09
1.169483637e-07
7.445202399e-07
3.876016309e-06
6.204312046e-07
5.089616453e-10
2.563342578e-13
-4.108963768e-09
+4.109081761499578e-09
8.534723944e-11
3.062921512e-14
5.084635637e-08
3.964438906e-09
9.56064042e-09
5.713458485e-07
8.602373058e-10
9.674381801e-09
5.575533604e-08
7.736937052e-10
5.25696461e-10
2.604611015e-08
5.185877553e-08
5.136914596e-07
4.457967612e-09
7.087532094e-06
2.866869842e-09
2.555695269e-07
1.14983411e-06
2.121499808e-11
6.538598731e-11
4.397576097e-10
1.68339254e-10
7.632177773e-10
3.194450213e-07
4.743002807e-08
5.174172592e-09
1.378394044e-06
0.0001075134541
0.0002104593459
1.291514926e-09
1.113812317e-09
2.836821935e-07
2.119414007e-05
3.257830135e-07
6.904331208e-08
1.133565101e-07
2.515736584e-07
2.868370103e-07
1.316255824e-07
1.143196584e-07
6.920318436e-10
1.261652011e-08
7.917924522e-13
7.954494796e-10
6.551135771e-09
5.085806996e-09
1.197408695e-11
7.765435603e-11
1.137649191e-10
9.984771335e-11
1.206636661e-07
9.837079112e-09
3.488285266e-08
4.035219127e-08
1.411962856e-08
1.030986206e-07
0.0001339305599
2.858760077e-12
3.300471596e-09
1.938419056e-10
3.176459434e-12
1.739161428e-09
4.024412744e-07
3.711166771e-08
3.68209028e-09
9.299354827e-11
4.060060302e-10
1.4612779e-05
9.474942894e-11
6.71405306e-14
2.422793691e-09
1.239892236e-07
2.586764905e-08
3.834829456e-09
4.981089193e-08
5.890823115e-10
1.712560509e-13
5.885874705e-12
3.850099718e-12
5.815334596e-10
2.871382336e-07
3.360539652e-09
4.474272729e-08
4.854260426e-09
6.835712438e-09
4.513662362e-10
1.937015923e-09
3.204387319e-07
2.834110865e-06
2.373908594e-10
8.959130596e-10
2.230101714e-07
2.935469223e-08
1.283653118e-05
1.768922904e-09
1.023937016e-07
1.810088219e-10
1.898013563e-06
1.555535859e-07
2.081843632e-11
1.130189635e-07
5.044203176e-10
4.419441096e-10
8.383860763e-08
1.081590063e-05
4.255989963e-09
4.961325451e-09
4.722044461e-12
2.605827829e-07
6.302277862e-12
7.401607027e-11
5.826773899e-09
3.775837499e-06
7.970709192e-06
1.077514227e-08
4.597223421e-10
6.830064994e-09
2.132411988e-10
3.596324757e-11
6.420547725e-09
2.577236401e-09
3.303229315e-09
6.18774819e-09
4.85351932e-08
4.830745338e-10
3.720076558e-06
1.356666145e-08
2.380000545e-09
6.934200897e-11
3.297569231e-09
3.435257108e-10
5.29774949e-08
4.382596642e-09
2.15867853e-08
5.682155814e-12
2.87290739e-07
1.68447438e-10
3.702374348e-08
4.342623412e-10
6.279260244e-08
2.403676166e-08
6.186547757e-09
1.906369995e-07
4.632234292e-10
2.755342608e-09
1.375328487e-08
2.872876599e-08
6.101397141e-07
3.268512848e-10
2.065059787e-07
2.629402879e-09
4.962476519e-11
1.388892596e-10
6.600566104e-08
2.769850168e-07
8.89254745e-07
4.790857584e-07
8.513434203e-09
2.371056219e-11
1.969882054e-08
1.194180775e-06
2.536127752e-06
3.538123306e-11
7.535922813e-11
3.279643366e-09
5.915578256e-09
3.02572337e-12
2.425014125e-08
1.736441207e-08
9.148487522e-11
2.773967859e-09
6.97620218e-08
1.281934525e-09
3.47740087e-11
7.026364385e-12
3.159133311e-08
2.31268605e-06
1.021090566e-11
6.311611218e-10
-3.496581068e-09
+3.496527402152061e-09
5.006407047e-09
1.086254132e-09
6.442837702e-07
3.40398573e-08
1.040533844e-07
2.015250069e-08
2.367399296e-08
1.455786943e-08
4.519557691e-07
6.28474478e-08
1.681943257e-08
3.350671592e-10
1.788948112e-11
8.659128848e-08
1.581729673e-10
3.16899147e-08
1.269452476e-09
3.307661094e-11
6.743047887e-11
4.30353612e-10
1.049230002e-07
7.726948943e-12
5.740327635e-07
7.011998453e-07
5.337360544e-09
6.364224996e-09
6.365307921e-10
7.724490255e-10
2.368415063e-08
-1.761682817e-07
+1.761719259305042e-07
3.074395438e-08
5.988468898e-06
1.432690031e-12
7.971925347e-11
1.851944696e-07
1.679623381e-10
6.480885684e-06
9.759515118e-09
3.448768806e-09
7.000972121e-11
2.170432861e-09
5.016191785e-09
1.49812676e-09
5.87915813e-10
6.048689475e-10
1.149696962e-08
8.392285319e-11
2.869985172e-13
9.120207147e-09
1.023039271e-06
1.843349283e-11
1.00493373e-06
1.437966503e-09
3.104285805e-06
3.736848159e-10
9.594801088e-11
1.232024101e-11
3.535931521e-07
7.131216184e-08
2.337326742e-10
2.82410779e-07
3.138271299e-08
5.494774408e-12
3.07020008e-12
2.011245542e-10
1.46145366e-11
3.082048934e-09
1.460968738e-10
3.186600527e-07
8.613894657e-10
4.227362638e-07
1.132945088e-08
1.385877833e-08
2.535124307e-09
5.647614453e-08
1.712204167e-10
9.689130473e-07
2.23306088e-09
2.607948681e-11
1.678634611e-13
2.258659145e-06
1.695727868e-09
1.706090682e-11
4.758749529e-12
1.426536131e-05
1.449437791e-07
4.244764945e-08
3.617494916e-10
8.345680248e-12
7.09135827e-13
8.179812515e-09
3.323888637e-09
1.436305431e-08
5.885631636e-10
3.584884074e-08
9.15498181e-09
9.640861937e-11
1.908907816e-09
2.558536978e-10
5.232929043e-10
1.261308954e-09
2.535469757e-10
-1.402765641e-06
+1.402741054722527e-06
9.405064753e-10
1.920855926e-07
4.40419547e-11
5.00514635e-10
1.406592654e-07
4.52773411e-10
4.146849421e-07
1.166300995e-07
7.051208905e-09
4.11857413e-09
9.269115967e-07
1.08066004e-10
1.045618377e-06
4.94758944e-07
6.160229109e-12
-1.351633908e-08
+1.351649588024934e-08
2.784256537e-09
1.934517844e-08
5.314666753e-13
1.357165448e-09
4.977436192e-12
6.477152561e-11
4.234578629e-08
4.918295292e-13
3.683451803e-08
3.727223271e-09
1.040244287e-09
6.604553745e-10
3.806952288e-09
1.097089218e-08
1.705682648e-08
6.319891143e-07
1.574130537e-09
3.013782963e-08
4.772880086e-09
3.679484278e-09
8.938611559e-10
5.601588127e-09
1.526952592e-05
3.497239887e-10
9.33354413e-08
1.359301132e-08
6.404370059e-11
1.535628496e-11
2.243494388e-09
4.662210325e-09
1.022052539e-10
1.383917366e-07
1.619072147e-12
1.507995029e-08
2.056138551e-09
1.984766447e-05
5.313284134e-09
1.04749306e-08
1.156452533e-07
1.201443255e-06
1.584111891e-09
2.09858386e-11
3.445445203e-09
4.164676e-06
4.092776954e-12
9.944337818e-10
5.276009823e-09
3.628978704e-11
2.197020798e-10
1.052815682e-08
1.068627029e-08
3.831431452e-11
1.364029236e-07
3.403124343e-10
5.38162347e-08
6.729323819e-12
1.207769187e-07
4.017542328e-09
1.315683749e-08
7.126836658e-11
1.491126652e-07
7.209030864e-10
6.628062787e-08
1.21591155e-10
6.230923352e-09
0.000101234359
5.30942773e-10
5.741385346e-10
2.429893426e-08
7.325998025e-12
3.533987787e-09
7.926231814e-10
8.291189826e-11
3.746302962e-11
5.242302307e-11
8.718447375e-08
6.176050515e-11
6.09961364e-09
2.644343829e-09
2.306901235e-08
7.322344218e-10
1.884045094e-06
3.755632286e-08
1.659778169e-08
5.495025114e-08
3.546727594e-09
2.674246643e-07
1.028252112e-09
9.289249703e-06
7.959960572e-08
1.068081079e-11
5.674343265e-07
6.214463157e-09
1.601066932e-06
1.977738539e-06
3.462786645e-07
7.373697384e-09
5.726726673e-12
diff --git a/t/ME_data/ME_Wp.dat b/t/ME_data/ME_Wp.dat
index 495c3be..7ce09cb 100644
--- a/t/ME_data/ME_Wp.dat
+++ b/t/ME_data/ME_Wp.dat
@@ -1,1289 +1,1289 @@
-0.0001355688725
+0.0001355719842764306
3.608722302e-05
4.306341977e-07
9.75550325e-07
7.965530581e-08
1.158741173e-06
2.661843554e-06
-2.305034582e-09
+2.30482384891786e-09
2.352135748e-07
1.82973514e-06
7.663338515e-10
1.356607225e-07
3.025163564e-07
9.674844388e-08
1.08577046e-08
1.140610128e-07
9.611408862e-06
9.926076721e-06
0.0002395617463
3.935060204e-07
1.135288401e-06
3.820198516e-06
9.707598452e-05
1.989871939e-08
4.670964333e-06
6.369460204e-08
5.769260099e-06
6.100771974e-07
4.177606056e-07
5.131150168e-05
2.623663346e-08
9.921267961e-07
2.81862487e-05
2.357728363e-06
1.52940914e-07
3.812119407e-09
4.089890436e-06
5.519579034e-08
3.108467376e-05
6.275083181e-07
0.0004711631212
1.064845939e-10
3.599907229e-07
0.0005693438579
5.344259311e-06
1.369519561e-07
4.610580273e-08
2.759136384e-08
0.0002011663129
8.79267475e-06
2.70389756e-06
1.871332754e-09
6.64520168e-07
0.0134591178
-8.026976255e-06
+8.025796760533041e-06
1.13118716e-07
2.178470183e-07
3.510890511e-05
1.899844473e-07
1.045126454e-06
1.869095429e-06
2.768546174e-06
8.551014192e-08
7.773400333e-07
4.53144298e-06
-1.601156732e-06
+1.601178083218329e-06
1.364165622e-07
5.96828212e-08
1.972891892e-07
7.870126436e-08
1.027083098e-09
1.68746816e-07
6.312089294e-10
3.391835464e-07
4.91589775e-06
1.471309784e-08
5.806965344e-07
1.107356095e-06
2.913951197e-05
2.124037647e-06
2.12990617e-07
1.356428765e-06
2.199225656e-08
9.711357262e-08
0.0001385574064
6.515223743e-07
2.535298939e-05
9.498820108e-06
7.275424133e-08
3.660521443e-08
1.057452388e-07
6.037059083e-08
1.535388646e-05
8.248395671e-07
3.681439971e-07
6.978801279e-07
0.0001134623229
0.00317235351
0.0002184140294
1.851577146e-06
2.387181833e-07
0.0002408737015
2.563590676e-07
1.00063893e-09
4.85334472e-06
1.03755407e-07
7.119588285e-06
2.561890817e-06
6.552784786e-08
2.501651241e-08
3.220116707e-08
4.669231392e-07
1.713317361e-05
1.245254933e-09
7.470170779e-06
1.12179479e-10
2.566806243e-08
0.000690252632
1.619690789e-08
0.007842918829
3.718578613e-05
0.0007152644659
3.366872929e-05
3.590280598e-05
0.0001083353837
1.746053028e-07
2.714474998e-05
5.49513934e-05
6.535082158e-05
5.041690215e-07
3.546092261e-06
1.20179547e-08
2.615186068e-08
0.005544506156
2.669153415e-08
1.124413241e-07
2.097675904e-07
0.0004497021745
3.522902236e-06
3.609795144e-10
2.381543307e-06
0.0001723836125
8.005109609e-06
2.55897517e-05
1.661304559e-05
7.830076033e-06
8.689122826e-05
1.60422181e-05
1.56327932e-07
7.475227738e-08
3.375304186e-08
7.208470216e-08
5.524153757e-10
5.4103721e-09
9.793444055e-05
1.21018852e-06
5.403237195e-06
9.396037827e-07
9.218690916e-07
1.916819566e-05
2.177871221e-06
0.0003020897612
6.197213888e-08
7.288410467e-08
1.364194249e-05
1.038835664e-05
1.181892855e-07
5.005165507e-06
3.611225632e-08
1.221974589e-08
1.942760585e-05
3.734309313e-06
0.0001093995311
2.068051727e-06
8.857501175e-06
4.173118322e-06
6.957949482e-09
1.44827703e-06
4.552053973e-07
1.371432834e-06
8.524886363e-07
1.659141663e-09
1.245863146e-07
0.0002242637793
3.699075972e-06
2.084644637e-10
2.212848552e-08
8.425415515e-08
2.259085562e-07
0.007259913649
4.467603509e-06
3.025351151e-06
1.224940717e-08
1.15612173e-07
1.910495778e-06
6.957098183e-08
3.814134721e-08
4.664811457e-07
1.051815649e-05
8.361499159e-09
6.131289047e-07
6.896616304e-06
2.90239901e-05
7.595782532e-09
1.528995487e-08
5.225540705e-09
3.429458529e-05
2.622252694e-09
0.0005682063516
7.453799242e-07
1.807444313e-05
9.740908615e-07
2.028479198e-06
9.039653406e-08
2.480442463e-06
0.0004804561921
1.082404311e-08
2.792579136e-06
0.0002121921202
9.812552393e-05
7.828388801e-07
0.001456303522
1.943547444e-08
8.760574342e-08
6.867301069e-10
1.796235147e-07
3.651695186e-07
1.227417157e-07
6.64659658e-06
8.406970187e-05
0.001885290476
0.0001570330885
2.314708659e-06
0.0002982727873
0.003229609179
2.42756184e-10
1.896039184e-06
1.190944746e-08
4.219946428e-07
7.272059576e-08
1.591962088e-06
2.437096009e-07
2.465458718e-10
3.867407921e-06
7.355333709e-09
3.882793753e-06
0.003249948478
1.636247528e-05
0.0003223792165
2.270764742e-05
1.158210386e-07
4.803440794e-06
1.938501144e-06
1.00011733e-06
7.966528549e-08
1.977164344e-05
3.697861147e-06
0.0001118273732
8.542125574e-09
3.549017984e-05
4.419576384e-08
0.000415723574
1.269135249e-06
1.336086492e-09
0.0001810768552
5.607526182e-08
7.073836906e-08
6.226377541e-05
3.842918719e-07
6.352284838e-05
1.767534549e-08
0.002554251878
1.646540099e-08
1.960645506e-08
0.0002646807702
1.199915361e-05
2.978511607e-08
0.000116967724
1.112499069e-08
4.229713001e-07
4.907009816e-08
1.380172172e-06
1.086563506e-05
4.639365797e-08
5.978111043e-07
7.468593061e-08
2.128107648e-08
3.477529199e-06
1.4537143e-09
3.088717592e-05
3.041172443e-05
2.604197881e-07
1.592445808e-06
1.544432896e-07
3.372379695e-06
0.0001468102323
5.574768414e-07
4.759127919e-10
0.00093067198
5.974247937e-06
0.0002033151473
2.248619959e-06
1.737481127e-05
3.508166199e-06
2.435498983e-09
6.369004222e-06
9.703729079e-09
1.740708919e-06
3.283306645e-07
1.701425138e-08
4.192916377e-08
4.498845691e-05
-3.834128247e-05
+3.834045529627708e-05
1.152621968e-09
1.27759416e-05
0.000281123799
1.725068848e-07
1.077501014e-06
0.004189517285
1.659706255e-07
1.343358849e-07
8.327570633e-07
2.346067359e-05
5.667318085e-07
1.495443328e-09
6.283129102e-07
7.691168852e-09
1.613371844e-08
4.367733457e-08
1.328415042e-09
1.452913836e-06
2.025565342e-07
1.29314441e-06
3.674495466e-07
-4.037738642e-06
+4.037694296051743e-06
2.086123437e-08
1.579986535e-06
7.231915707e-10
3.157015643e-07
0.0003660037209
0.0008949054777
5.638174345e-08
3.288585027e-06
5.432504634e-08
5.180759332e-10
5.121204016e-08
1.939451095e-08
1.979142335e-05
2.529002486e-08
0.0002538283082
3.80397913e-05
7.172326439e-06
2.759121311e-06
3.708767178e-06
4.71298987e-08
0.0001084418497
5.277822142e-08
4.364063475e-08
1.209897825e-07
3.136890863e-08
1.044095011e-07
6.962956326e-08
6.454347396e-06
2.03105437e-09
4.064764112e-07
1.125205343e-06
3.916697214e-06
5.611785349e-07
1.848611741e-07
2.582403728e-07
8.839788397e-05
2.850326226e-06
2.108292593e-07
1.232805123e-08
9.539264139e-08
0.006539873531
-2.250970922e-05
+2.250946039109731e-05
1.89348364e-07
2.619177604e-11
3.531432382e-10
1.155035491e-07
1.892464014e-07
1.339909548e-08
5.52939388e-08
1.151900346e-05
0.005889722048
8.216256879e-05
1.557120111e-07
9.089312296e-08
4.410189849e-07
4.12565868e-08
1.888812372e-07
5.726907174e-05
1.874558888e-09
9.008463668e-12
1.709044906e-05
1.67931414e-07
8.335566347e-06
1.14159826e-08
0.002680456252
1.068398085e-08
1.957173741e-06
2.244367604e-05
1.231910962e-05
4.629935916e-05
4.576226372e-06
1.980987784e-08
3.93911073e-07
2.650096111e-05
8.453749094e-06
3.871976822e-08
0.002862691182
4.135291611e-08
0.0002240310514
1.626024307e-07
1.446697749e-07
1.81572275e-06
3.82343645e-07
1.519921801e-06
7.180741242e-10
3.48543689e-05
7.358200622e-08
7.791803572e-06
7.340099307e-07
1.64882977e-05
5.701378479e-07
1.101811989e-07
2.180113388e-08
0.005048614688
1.4109536e-08
1.796629531e-05
3.537585755e-05
2.035011121e-05
1.003894855e-08
0.0007063054044
8.670678472e-05
6.82639369e-08
5.900706666e-06
1.447352541e-08
1.072050501e-05
2.460311e-05
3.389946454e-08
3.777468751e-07
7.004890999e-07
3.633307241e-05
4.394910201e-07
7.401664237e-07
0.004667703381
4.822543865e-07
1.289604844e-08
0.0001296635825
1.387032601e-06
2.915976119e-05
1.716214639e-07
2.901236713e-08
7.505249455e-10
7.295054917e-08
5.765824783e-07
3.265782722e-08
1.700596863e-06
3.410463387e-08
7.817851611e-06
0.000171359498
4.092934015e-08
9.55513037e-08
2.412195475e-06
3.183516928e-09
0.0006930495435
1.228003613e-07
8.804634278e-05
2.13517622e-06
4.693017624e-06
5.984372321e-08
2.487744823e-06
1.776206816e-07
6.533684691e-07
1.608251404e-06
4.483048832e-06
1.390213091e-06
2.463272722e-09
0.0001373179166
4.692519624e-10
6.223649822e-09
3.026899135e-06
5.540599649e-08
5.422865018e-08
8.715360663e-08
3.005262078e-06
4.491518597e-06
1.296331911e-08
8.034189259e-07
8.445679444e-07
1.331376027e-08
-1.808292903e-06
+1.808252927610846e-06
4.811816349e-06
3.356636075e-06
3.718204067e-07
1.63321067e-06
1.145748524e-06
1.082083666e-06
6.645503138e-07
7.690061877e-06
9.446655788e-07
0.0001628166399
0.0004161528763
2.733182188e-06
4.797973294e-06
0.0001243192011
0.0004489774904
4.018994469e-08
6.294131882e-05
0.0006300938801
3.276964999e-07
5.253447013e-07
6.013623386e-09
6.133379962e-07
7.498082273e-08
9.70785668e-08
2.892020441e-07
0.0008858500922
5.834224897e-12
2.657258192e-05
1.937050603e-08
1.014207639e-08
2.078416478e-09
9.09266244e-09
1.171070677e-06
6.186008954e-07
7.316626453e-08
1.137928043e-05
3.926244745e-07
9.411347464e-08
2.358883261e-05
4.562416832e-07
8.746279414e-07
2.158976175e-07
7.682516501e-08
4.623476813e-08
1.822643012e-09
5.649636691e-07
8.861102166e-07
3.661799027e-06
0.005151284527
1.333261465e-08
0.01070078066
3.557312307e-07
7.956062001e-07
1.102026215e-05
9.099085434e-08
9.760469392e-07
4.026083515e-07
3.148718342e-05
3.652137087e-08
5.83259011e-08
0.000275860311
5.979810844e-06
2.336312158e-06
1.214278297e-06
1.257019376e-05
4.380094449e-07
8.01469852e-07
0.003169580958
9.221636814e-05
4.690760304e-07
7.648050266e-07
1.418690988e-05
2.220802016e-07
3.14201845e-07
1.911591064e-07
2.452984351e-08
0.0001445594492
2.06149244e-07
1.315444637e-09
3.270418535e-05
2.302140766e-07
4.002096106e-06
1.609778389e-06
1.156198845e-08
4.49026944e-08
0.0002105189957
1.502882545e-10
0.0003938508954
2.359522705e-06
5.766519564e-05
4.012744154e-07
1.060136498e-07
0.0002141768308
2.679472087e-06
1.069266147e-05
2.48432308e-09
0.0001019831675
1.193130845e-06
1.086711459e-06
1.985358343e-07
1.172768569e-10
4.982454227e-07
2.574857412e-05
0.001250253279
2.058059351e-09
1.333980604e-07
9.818799108e-09
5.15408822e-07
0.0001083932826
1.428558882e-06
1.74625784e-05
2.151559189e-05
0.01028504249
7.324000358e-07
2.685867501e-08
4.250514247e-05
0.0004356795916
3.289457433e-07
0.00323686003
3.007561717e-08
2.793887221e-06
1.897911806e-05
1.61972882e-07
8.953729168e-05
4.802035152e-06
7.858920744e-06
7.610825707e-07
2.314286712e-05
4.93977485e-09
4.92614674e-09
8.469370311e-07
3.80925189e-05
8.402086849e-09
4.426792197e-07
2.426545735e-06
3.814110561e-08
1.746240227e-06
1.52446212e-07
7.921180881e-05
2.10135866e-07
4.839236788e-07
5.352097408e-07
4.622274946e-05
0.001900140645
1.924092071e-06
0.001430541179
4.980473899e-08
0.001646505316
1.600615719e-05
2.692315114e-06
4.189143283e-06
9.988552464e-08
3.325410404e-05
1.063381388e-06
4.230272945e-06
3.9885381e-06
0.0001050685482
5.986074981e-08
2.025564094e-06
4.345191943e-07
8.910595474e-05
2.543279105e-06
0.0002129944726
8.738486781e-06
4.913962453e-09
3.3533254e-08
3.330632607e-05
1.95151118e-06
1.384784661e-05
2.148337293e-07
1.685318586e-07
1.114845673e-07
6.020836887e-05
1.280889433e-07
4.215297019e-07
3.062431208e-07
6.986710594e-07
1.615334662e-07
1.317663183e-07
0.0008030217322
2.015927418e-05
5.19186806e-06
5.500044264e-07
7.15420681e-08
0.0002157285841
9.440044669e-07
6.51295213e-08
1.03903577e-08
0.0001795918463
1.962893454e-05
1.739516144e-07
2.12576483e-08
1.743490745e-05
1.541416437e-06
3.391980216e-06
8.359096999e-06
1.415613908e-09
7.887561636e-07
1.745081946e-07
1.932727685e-07
9.828062711e-07
1.291461662e-07
9.278171187e-05
4.264104098e-06
3.754457748e-07
3.452634763e-05
8.101309385e-06
4.056256865e-09
2.25304955e-07
1.172523704e-05
4.796449304e-07
3.785295788e-09
2.390481542e-05
5.308836572e-08
4.670595432e-05
6.758171771e-05
3.147635808e-09
1.661568152e-05
5.893560303e-08
3.900939999e-05
3.504262189e-08
8.675964674e-05
0.0001335992511
1.302314855e-06
1.146394419e-06
6.394209116e-06
2.519496842e-05
0.0001367050958
4.154554201e-09
2.099823438e-06
2.146448613e-05
1.233863159e-09
-0.001728586801
+0.001728401026344979
4.09215386e-05
0.007151464488
4.760189124e-09
0.0007525692076
1.763213602e-08
2.706320225e-07
4.168939514e-06
6.813774053e-07
6.199680386e-06
1.450358601e-09
0.001178944008
3.552524877e-07
0.001006338796
3.68488544e-07
2.167674892e-07
2.537842791e-06
7.658583099e-07
2.591411922e-07
4.604369581e-08
4.168373998e-07
7.238778845e-09
1.170127806e-06
4.393582652e-07
0.0007919330066
1.02877764e-06
6.918341723e-08
9.436339417e-08
1.591183988e-06
2.089478745e-06
1.613369096e-07
3.637102293e-07
9.604977934e-07
9.357975878e-07
2.88751538e-06
6.126247948e-06
1.015660592e-05
5.755046456e-06
0.002323684387
8.012994582e-07
5.26519905e-07
0.0002032569161
1.451031711e-07
5.477022502e-07
7.468911459e-07
1.172272701e-05
5.000191366e-06
1.572889827e-06
6.442064306e-08
2.908881721e-05
4.569175874e-05
0.0003671971117
1.009483851e-07
1.696662855e-09
1.138346461e-08
7.716922378e-07
3.529252959e-06
6.363107685e-07
1.932801483e-07
-5.985897033e-07
+5.98583507415228e-07
3.856271125e-06
8.483829284e-07
7.519170568e-09
6.046535116e-06
8.332016051e-05
-0.0004929003218
+0.0004929270085709528
1.066348909e-06
6.354668796e-05
2.325179903e-05
3.330122089e-07
7.907070499e-09
3.514991564e-06
6.69519932e-06
7.791955281e-08
3.007406815e-09
1.687262337e-06
2.966015995e-06
2.551884972e-08
-0.02757085539
+0.02757047144357582
1.902124521e-06
1.858863499e-06
1.174589642e-07
6.265121517e-07
0.0001498056382
2.510077007e-07
0.003704558135
7.86457155e-05
1.776365133e-08
0.0002011705313
4.239853347e-06
1.056573955e-06
0.0008521527377
4.768494505e-07
5.244547978e-06
-0.0001425058496
+0.0001425035783425772
4.507288271e-07
3.092538831e-07
0.007681533942
1.571813043e-07
3.341399818e-08
9.577252789e-07
1.135356146e-09
1.20086557e-07
4.167792379e-06
7.827764185e-06
1.444125375e-06
2.567040739e-05
3.703109657e-07
2.683391447e-05
2.223174582e-07
1.771598439e-05
1.154475457e-07
4.676107532e-07
1.045592241e-08
5.174223151e-07
8.734482888e-07
1.146406702e-08
8.784052908e-07
7.500525642e-10
1.921656528e-05
5.899142263e-08
7.30476082e-05
1.070817446e-08
2.432332991e-05
2.466303926e-09
1.587801136e-05
5.253580132e-08
1.186388137e-06
2.513054181e-08
2.878235795e-05
1.268594957e-09
4.036344626e-07
1.480245097e-05
0.007764649525
1.266626553e-05
2.673316923e-07
2.052285557e-07
3.026411522e-09
3.958784188e-09
5.153150957e-10
2.304903761e-10
4.00173635e-09
5.422672396e-12
1.0797419e-08
7.876103608e-07
2.04929605e-13
2.452584831e-10
6.723994081e-09
1.881390676e-09
8.769218159e-09
1.325727129e-07
8.79685454e-10
8.598804962e-09
1.838865267e-09
1.411793749e-09
3.180735014e-09
1.657163946e-09
6.291662125e-07
1.407385591e-08
2.446407045e-10
3.694683988e-06
7.168303932e-08
1.515131945e-07
5.976833117e-07
9.943757029e-11
1.770333656e-10
6.826832537e-10
3.501999566e-09
1.003726202e-11
-2.246655664e-10
+2.246679297996575e-10
1.56313715e-08
1.258803718e-08
1.510698504e-09
1.860286446e-08
7.436041453e-11
8.891989836e-11
8.343726115e-07
7.82113065e-10
7.372199042e-09
1.73794743e-08
4.796628283e-10
1.006843364e-11
2.297522292e-07
2.964339936e-09
5.119472399e-06
1.217098723e-09
1.247855962e-10
2.379923838e-10
6.847023572e-08
1.126084853e-08
9.36637247e-09
8.815353934e-13
1.248942677e-09
1.288006978e-06
1.862006179e-05
5.986024169e-09
3.576379278e-09
1.396207867e-07
5.597776909e-09
1.227869194e-07
3.184774841e-10
1.281546379e-09
3.785370242e-08
8.065384687e-10
1.797326593e-08
4.870914289e-09
-1.481397924e-08
+1.481421747250364e-08
6.148653826e-13
2.60768604e-08
7.503399778e-09
7.319978311e-10
3.917579972e-06
2.456000696e-09
1.6358728e-09
1.504286785e-07
4.626076109e-12
2.544705854e-10
2.326334049e-09
1.209466619e-07
9.766234969e-09
2.844168649e-07
1.201220077e-06
7.29513531e-08
2.673472495e-08
3.217664769e-11
3.051103391e-08
8.284753217e-08
2.058274741e-10
2.829620954e-09
6.31774405e-05
5.647654886e-08
3.939849241e-09
6.945357028e-09
2.688851327e-09
1.246699666e-07
1.659903866e-07
2.236198294e-05
9.480526263e-09
9.500853178e-07
7.271132888e-08
4.043174967e-07
7.378138213e-10
1.151223659e-07
3.256246906e-08
1.170385775e-12
5.503748995e-08
1.04382963e-11
1.203555924e-08
1.075568292e-10
8.108085296e-07
4.847854101e-11
7.343632975e-06
1.187566413e-11
3.504032934e-08
7.018111299e-12
2.463652527e-10
1.851835159e-10
1.513115041e-09
8.018679868e-09
3.929445577e-06
2.084075835e-07
5.951216671e-10
3.507041185e-06
5.267975954e-12
3.372609083e-06
1.694025855e-11
1.867026183e-10
2.116504705e-10
3.831483929e-10
6.969465117e-05
1.539468292e-08
1.101455262e-11
1.388260351e-10
5.247893014e-08
4.682122453e-08
1.641630928e-10
2.19360099e-08
1.553226092e-09
2.996160973e-08
2.723720456e-11
1.94231476e-09
1.292957006e-11
3.828519911e-12
1.610161407e-08
1.687340225e-06
1.925507425e-09
4.366763058e-13
6.454545293e-07
6.024146742e-07
2.580279916e-10
1.225824367e-05
1.999348472e-11
4.860969609e-11
2.362388708e-08
1.373681959e-08
3.790417277e-12
1.087174385e-08
3.1293099e-07
1.364398624e-09
1.715330241e-05
3.172978773e-08
1.823609884e-08
3.114529471e-09
3.261473002e-10
2.126205049e-05
5.67083409e-08
9.666274841e-10
7.465321603e-10
1.354621048e-10
2.460861072e-08
1.087512798e-09
2.736687033e-07
5.113035943e-10
6.306214627e-07
3.166204088e-09
2.044825034e-08
6.895612973e-10
6.672832434e-11
5.398174211e-09
1.015670208e-07
5.845735119e-10
4.836785182e-08
5.100824595e-10
1.572989395e-08
1.787382013e-10
1.632652152e-08
7.01728011e-08
7.320289671e-12
3.029928398e-11
1.409128876e-12
1.238566734e-06
8.502740668e-08
3.438552357e-08
4.863721359e-07
5.154442037e-07
2.08972413e-09
3.431301186e-10
1.781922807e-08
1.544285017e-09
4.30867667e-09
2.652046726e-09
-1.604339305e-11
+1.604320243519326e-11
1.019689627e-05
1.123301038e-09
2.815571339e-12
1.072295152e-10
6.773569516e-09
2.399435949e-09
1.754778047e-07
6.118111465e-11
8.592783136e-06
6.887731648e-11
1.799971692e-09
2.674451282e-10
1.267736949e-08
1.616784109e-09
2.41843279e-10
2.161775938e-07
3.594767539e-06
1.314807907e-05
3.035193063e-06
5.146951027e-12
8.12518692e-07
2.904383885e-07
8.278942839e-08
1.491259815e-12
3.232983294e-10
3.014861189e-07
3.869348398e-08
4.060145232e-09
2.722588887e-09
1.252225008e-06
8.534088341e-10
1.170269802e-13
1.013048217e-07
1.173680983e-10
1.230880253e-06
3.603057796e-07
2.340496864e-10
1.622753822e-07
1.855111476e-10
1.873444017e-13
1.477986621e-10
5.974804124e-09
5.502900798e-08
6.238803852e-09
2.281648706e-08
2.708668262e-09
8.942323868e-07
8.600240612e-08
5.313999062e-08
1.150328412e-06
4.275682161e-10
2.535797082e-12
1.472728159e-07
2.05454394e-10
2.452447725e-06
2.149816554e-09
1.302067004e-08
1.354519139e-05
1.493773984e-10
4.155192241e-09
1.975441615e-09
7.883829674e-09
3.017038977e-08
6.10674516e-07
-7.43296019e-09
+7.433047934441598e-09
3.802076048e-08
1.356542979e-08
7.945938899e-11
6.146929649e-09
1.061272859e-06
3.425627163e-09
4.724471432e-09
2.264988832e-06
5.206993128e-12
1.903839604e-11
5.585175902e-09
5.143392872e-12
8.796344869e-08
2.018357672e-09
7.316136952e-07
6.365052389e-11
5.084402628e-08
4.955636143e-10
7.579083396e-10
1.041589492e-10
1.067922679e-10
2.629804234e-10
7.253503328e-09
3.644406678e-09
3.524116213e-06
1.911509485e-10
-7.70936902e-09
+7.705863931467e-09
9.497529697e-06
1.776186622e-07
7.672088841e-12
1.053749536e-07
4.528511033e-07
2.689071764e-11
2.172737406e-09
1.388580447e-10
8.063951339e-10
5.664084936e-08
1.744805297e-09
2.052101897e-11
2.110860659e-11
4.471965163e-07
1.547872897e-09
1.123990826e-08
1.220466442e-08
7.869401405e-07
1.643506226e-10
1.008594225e-08
3.353990229e-07
6.498833539e-10
1.784315183e-10
3.805120267e-07
8.427310035e-07
1.498076756e-09
1.039474025e-10
1.111563184e-08
1.747513597e-09
2.355165317e-09
1.678233595e-08
2.202228399e-10
5.649816397e-08
7.451025295e-08
3.278979106e-11
3.654252804e-11
1.935635959e-11
4.454493836e-09
1.468843362e-11
6.571277705e-09
6.056988847e-10
1.451019442e-07
1.129031793e-10
1.042132383e-11
3.257388686e-11
9.125965424e-09
1.157787061e-08
1.347062065e-10
1.77255365e-07
2.066007246e-07
1.183381206e-10
4.34066107e-11
1.724089241e-07
1.017527409e-09
4.44819419e-08
9.298034227e-10
1.457982094e-08
1.344707431e-06
3.942998333e-10
6.431805784e-10
2.598344406e-09
9.723411144e-11
5.888094133e-07
1.250960845e-08
1.918965372e-08
1.266172515e-09
2.358561388e-08
9.312555038e-07
4.777092354e-09
1.101320525e-05
4.464641617e-10
7.607867556e-08
5.154678729e-10
2.264518558e-08
8.764203667e-09
8.615359173e-07
1.400305016e-06
2.971220423e-07
1.703068601e-09
9.085894816e-08
7.676840903e-09
3.495641886e-07
5.3486215e-10
8.550499106e-10
5.359621176e-10
3.270585775e-08
2.199341422e-11
3.481833841e-07
6.965999199e-10
1.079185990617177e-09
2.042388686e-11
1.522460286e-08
1.814687048e-09
3.79827064e-08
3.251582144e-06
1.339564566e-08
1.067207794e-10
2.907408741e-09
2.121056487e-09
1.929150342e-10
5.485146663e-10
2.944552971e-10
-3.400695511e-09
+3.400562267759813e-09
2.722315508e-07
1.103475015e-07
4.918180357e-10
1.403530754e-08
1.186111191e-10
2.383851903e-09
6.366433869e-10
2.884801368e-08
3.233653883e-09
0.0001662977498
2.949009865e-10
2.818141905e-08
2.146479379e-09
2.878567552e-11
1.429984095e-11
2.151490092e-11
4.991815907e-07
7.608561997e-11
4.322238014e-08
-8.515965863e-10
+8.515744986817629e-10
4.522028585e-09
4.658966814e-08
4.872474226e-10
1.225955173e-08
9.732456294e-09
5.55287071e-08