Commit 6e1c06e9 authored by Marek Schoenherr's avatar Marek Schoenherr
Browse files

Higgs core scale setter (and nothing else this time)

parent 82943992
#include "PHASIC++/Scales/Core_Scale_Setter.H"
#include "PHASIC++/Process/Single_Process.H"
#include "ATOOLS/Org/Exception.H"
#include "ATOOLS/Org/Message.H"
namespace PHASIC {
class Higgs_Core_Scale: public Core_Scale_Setter {
public:
Higgs_Core_Scale(const Core_Scale_Arguments &args):
Core_Scale_Setter(args) {}
PDF::CParam Calculate(ATOOLS::Cluster_Amplitude *const ampl);
ATOOLS::Cluster_Amplitude *Cluster
(ATOOLS::Cluster_Amplitude *const ampl) const;
};// end of class Scale_Setter_Base
}// end of namespace PHASIC
using namespace PHASIC;
using namespace ATOOLS;
PDF::CParam Higgs_Core_Scale::Calculate(Cluster_Amplitude *const ampl)
{
DEBUG_FUNC("");
msg_Debugging()<<*ampl<<"\n";
if (ampl->Legs().size()==3 && ampl->NIn()==2) {
double kt2cmin(ampl->Leg(2)->Mom().Abs2());
return PDF::CParam(kt2cmin,kt2cmin,0.0,kt2cmin,-1);
}
double muf2(0.0), mur2(0.0), muq2(0.0);
Cluster_Amplitude *campl(Cluster(ampl->Copy()));
if (campl->Legs().size()!=ampl->Legs().size())
msg_Debugging()<<*campl<<"\n";
if (campl->Legs().size()!=4) {
double q2((campl->Leg(0)->Mom()+campl->Leg(1)->Mom()).Abs2());
Vec4D ewsum;
for (size_t i(0);i<campl->Legs().size();++i)
if (!campl->Leg(i)->Flav().Strong()) ewsum+=campl->Leg(i)->Mom();
if (ewsum==Vec4D()) ewsum=campl->Leg(0)->Mom()+campl->Leg(1)->Mom();
if (campl->NIn()==2 &&
campl->Leg(0)->Flav().Strong() &&
campl->Leg(1)->Flav().Strong()) {// HThat'/2
q2=ewsum.PPerp();
for (size_t i(0);i<campl->Legs().size();++i)
if (campl->Leg(i)->Flav().Strong())
q2+=campl->Leg(i)->Mom().PPerp();
q2=sqr(ewsum.Mass()+q2/2.0);
}
campl->Delete();
return PDF::CParam(q2,dabs(ewsum.Abs2()),0.0,q2,-1);
}
Flavour fl[4]={campl->Leg(0)->Flav(),campl->Leg(1)->Flav(),
campl->Leg(2)->Flav(),campl->Leg(3)->Flav()};
if (fl[0].Strong() && fl[1].Strong()) {// hh collision
if (fl[2].Strong() && fl[3].Strong()) {
msg_Debugging()<<"pure QCD like\n";
double s(2.0*campl->Leg(0)->Mom()*campl->Leg(1)->Mom());
double t(2.0*campl->Leg(0)->Mom()*campl->Leg(2)->Mom());
double u(2.0*campl->Leg(0)->Mom()*campl->Leg(3)->Mom());
muq2=muf2=mur2=-1.0/(1.0/s+1.0/t+1.0/u);
}
else if (!fl[2].Strong() && !fl[3].Strong()) {
msg_Debugging()<<"DY like\n";
muq2=muf2=mur2=(campl->Leg(0)->Mom()+campl->Leg(1)->Mom()).Abs2();
}
else if (fl[2].Strong() && !fl[3].Strong()) {
msg_Debugging()<<"jV like\n";
muq2=muf2=mur2=campl->Leg(3)->Mom().MPerp2();
}
else if (!fl[2].Strong() && fl[3].Strong()) {
msg_Debugging()<<"Vj like\n";
muq2=muf2=mur2=campl->Leg(2)->Mom().MPerp2();
}
else THROW(not_implemented,"Please define your own core scale");
}
else if (!fl[0].Strong() && !fl[1].Strong()) {// ll collision
if (fl[2].Strong() && fl[3].Strong()) {
msg_Debugging()<<"jets like\n";
muq2=muf2=mur2=(campl->Leg(0)->Mom()+campl->Leg(1)->Mom()).Abs2();
}
}
else {// lh collision
THROW(not_implemented,"Please define your own core scale");
}
campl->Delete();
msg_Debugging()<<"\\mu_f = "<<sqrt(muf2)<<"\n"
<<"\\mu_r = "<<sqrt(mur2)<<"\n"
<<"\\mu_q = "<<sqrt(muq2)<<"\n";
return PDF::CParam(muf2,muq2,0.0,mur2,-1);
}
Cluster_Amplitude *Higgs_Core_Scale::Cluster
(Cluster_Amplitude *const ampl) const
{
if (ampl->Legs().size()==ampl->NIn()+2) return ampl;
Single_Process *proc(ampl->Proc<Single_Process>());
for (size_t i(ampl->NIn());i<ampl->Legs().size();++i) {
Cluster_Leg *li(ampl->Leg(i));
for (size_t j(i+1);j<ampl->Legs().size();++j) {
Cluster_Leg *lj(ampl->Leg(j));
if (proc->Combinable(li->Id(),lj->Id())) {
Flavour_Vector fls(proc->CombinedFlavour(li->Id()|lj->Id()));
for (size_t k(0);k<fls.size();++k) {
bool dec(false);
for (size_t l(0);l<ampl->Decays().size();++l)
if (ampl->Decays()[l]->m_id==(li->Id()|lj->Id())) {
dec=true;
break;
}
if ((!li->Flav().Strong() && !lj->Flav().Strong() &&
!fls[k].Strong()) || dec) {
msg_Debugging()<<"combine "<<ID(li->Id())<<"&"<<ID(lj->Id())
<<"->"<<fls[k]<<" ("<<dec<<")\n";
li->SetFlav(fls[k]);
li->SetMom(li->Mom()+lj->Mom());
li->SetId(li->Id()|lj->Id());
lj->Delete();
for (ClusterLeg_Vector::iterator lit(ampl->Legs().begin());
lit!=ampl->Legs().end();++lit)
if (*lit==lj) {
ampl->Legs().erase(lit);
break;
}
return Cluster(ampl);
}
}
}
}
}
return ampl;
}
DECLARE_ND_GETTER(Higgs_Core_Scale,"HIGGS",
Core_Scale_Setter,Core_Scale_Arguments,true);
Core_Scale_Setter *ATOOLS::Getter
<Core_Scale_Setter,Core_Scale_Arguments,Higgs_Core_Scale>::
operator()(const Core_Scale_Arguments &args) const
{
return new Higgs_Core_Scale(args);
}
void ATOOLS::Getter<Core_Scale_Setter,Core_Scale_Arguments,
Higgs_Core_Scale>::
PrintInfo(std::ostream &str,const size_t width) const
{
str<<"Higgs core scale";
}
......@@ -9,6 +9,7 @@ libPhasicScales_la_SOURCES = \
Scale_Setter_Base.C \
Core_Scale_Setter.C \
Default_Core_Scale.C \
Higgs_Core_Scale.C \
KFactor_Setter_Base.C \
Tag_Setter.C \
Variable_Scale_Setter.C \
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment