// Channel_Generator3V #include "PHASIC++/Channels/Single_Channel.H" #include "ATOOLS/Org/Run_Parameter.H" #include "ATOOLS/Org/MyStrStream.H" #include "ATOOLS/Org/Scoped_Settings.H" #include "PHASIC++/Channels/Channel_Elements.H" #include "PHASIC++/Channels/Vegas.H" using namespace PHASIC; using namespace ATOOLS; namespace PHASIC { class C3_0 : public Single_Channel { double m_thexp; double m_amct,m_alpha,m_ctmax,m_ctmin; Info_Key m_kTC23_0_4_1__2_3,m_kTC23_0__1__4_23,m_kZS_0; Vegas* p_vegas; public: C3_0(int,int,Flavour*,Integration_Info * const); ~C3_0(); void GenerateWeight(Vec4D *,Cut_Data *); void GeneratePoint(Vec4D *,Cut_Data *,double *); void AddPoint(double); void MPISync() { p_vegas->MPISync(); } void Optimize() { p_vegas->Optimize(); } void EndOptimize() { p_vegas->EndOptimize(); } void WriteOut(std::string pId) { p_vegas->WriteOut(pId); } void ReadIn(std::string pId) { p_vegas->ReadIn(pId); } void ISRInfo(int &,double &,double &); std::string ChID(); }; } extern "C" Single_Channel * Getter_C3_0(int nin,int nout,Flavour* fl,Integration_Info * const info) { return new C3_0(nin,nout,fl,info); } void C3_0::GeneratePoint(Vec4D * p,Cut_Data * cuts,double * _ran) { double *ran = p_vegas->GeneratePoint(_ran); for(int i=0;iGetscutAmegic(std::string("23")); Vec4D p23; double s23 = CE.ThresholdMomenta(m_thexp,4.*sqrt(s23_min),s23_min,s23_max,ran[0]); double s4 = ms[4]; double tmass201 = Flavour((kf_code)(23)).Mass(); m_ctmax = cuts->cosmax[0][4]; m_ctmin = cuts->cosmin[0][4]; CE.TChannelMomenta(p[0],p[1],p[4],p23,s4,s23,tmass201,m_alpha,m_ctmax,m_ctmin,m_amct,0,ran[1],ran[2]); double tmass202 = Flavour((kf_code)(23)).Mass(); Vec4D p0_4 = p[0]-p[4]; CE.TChannelMomenta(p0_4,p[1],p[2],p[3],s2,s3,tmass202,m_alpha,1.,-1.,m_amct,0,ran[3],ran[4]); } void C3_0::GenerateWeight(Vec4D* p,Cut_Data * cuts) { double wt = 1.; Vec4D p234=p[0]+p[1]; double s234_max = p234.Abs2(); double s23_max = sqr(sqrt(s234_max)-sqrt(ms[4])); double s3 = ms[3]; double s2 = ms[2]; double s23_min = cuts->GetscutAmegic(std::string("23")); Vec4D p23 = p[2]+p[3]; double s23 = dabs(p23.Abs2()); wt *= CE.ThresholdWeight(m_thexp,4.*sqrt(s23_min),s23_min,s23_max,s23,rans[0]); double s4 = ms[4]; double tmass201 = Flavour((kf_code)(23)).Mass(); m_ctmax = cuts->cosmax[0][4]; m_ctmin = cuts->cosmin[0][4]; if (m_kTC23_0__1__4_23.Weight()==ATOOLS::UNDEFINED_WEIGHT) m_kTC23_0__1__4_23<GenerateWeight(rans); if (wt!=0.) wt = vw/wt/pow(2.*M_PI,3*3.-4.); weight = wt; } C3_0::C3_0(int nin,int nout,Flavour* fl,Integration_Info * const info) : Single_Channel(nin,nout,fl) { Settings& s = Settings::GetMainSettings(); name = std::string("C3_0"); rannum = 5; rans = new double[rannum]; m_thexp = s["THRESHOLD_EPSILON"].Get(); m_amct = 1.0 + s["CHANNEL_EPSILON"].Get(); m_alpha = s["TCHANNEL_ALPHA"].Get(); m_ctmax = 1.; m_ctmin = -1.; m_kTC23_0_4_1__2_3.Assign(std::string("TC23_0_4_1__2_3"),2,0,info); m_kTC23_0__1__4_23.Assign(std::string("TC23_0__1__4_23"),2,0,info); m_kZS_0.Assign(std::string("ZS_0"),2,0,info); p_vegas = new Vegas(rannum,100,name); } C3_0::~C3_0() { delete p_vegas; } void C3_0::ISRInfo(int & type,double & mass,double & width) { type = 2; mass = 0; width = 0.; } void C3_0::AddPoint(double Value) { Single_Channel::AddPoint(Value); p_vegas->AddPoint(Value,rans); } std::string C3_0::ChID() { return std::string("CG2$MTH_23$TC23_0_4_1__2_3$TC23_0__1__4_23$ZS_0$"); }