#include "ATOOLS/Phys/Flavour.H" #include "HADRONS++/PS_Library/ResonanceFlavour.H" #include "ATOOLS/Org/Message.H" #include "ATOOLS/Org/My_File.H" #include "ATOOLS/Org/My_MPI.H" #include "ATOOLS/Org/Shell_Tools.H" #include "ATOOLS/Org/Run_Parameter.H" using namespace HADRONS; using namespace ATOOLS; using namespace std; // 3 particle phase space function lambda double ResonanceFlavour::Lambda( double a, double b, double c ) { double L = sqr(a-b-c)-4.*b*c; if (L>0.) return L; return 0.; } double ResonanceFlavour::Sqrt_Lambda( double a, double b, double c) { return sqrt(Lambda(a,b,c)); } SimpleResonanceFlavour::SimpleResonanceFlavour( std::string name, double _mass, double _width ) : m_name(name), m_mass(_mass), m_width(_width), m_mass2(sqr(_mass)) { } ResonanceFlavour::ResonanceFlavour( kf_code _kfc, double _mass, double _width, int _run ) : SimpleResonanceFlavour( Flavour(_kfc).IDName(), _mass, _width) { m_kfc = _kfc; m_running = _run; m_body = 2; p_hist = NULL; m_G_at_m2 = 1.; } void ResonanceFlavour::InitialiseThreeBodyResonance( ResonanceFlavour &res1 ) { InitialiseThreeBodyResonance(res1,res1,0.); } void ResonanceFlavour::InitialiseThreeBodyResonance( ResonanceFlavour &res1, ResonanceFlavour &res2, double beta ) { m_body=3; if( m_running ) { p_hist = CreateGHistogram( res1, res2, beta, kf_pi_plus ); m_G_at_m2 = GetValueOfG(m_mass2); } } Histogram * ResonanceFlavour::CreateGHistogram( ResonanceFlavour res1, ResonanceFlavour res2, double beta, kf_code out ) { // create file name char fn[512]; snprintf(fn, 512, "GQ2_Mres=%.3f_Gres=%.3f_MresP=%.3f_GresP=%.3f_beta=%.3f_Mout=%.3f.dat", res1.Mass(), res1.Width(), res2.Mass(), res2.Width(), beta, Flavour(out).HadMass() ); // look if file already exists std::string filename=rpa->gen.Variable("SHERPA_SHARE_PATH")+"/PhaseSpaceFunctions/"+string(fn); // if file does not exist if (!FileExists(filename)) { // create histogram (i.e. table of values) msg_Out()<<"Create necessary phase space function for chosen parameters.\n" <<"This may take some time. Please wait..."<Insert( q2, phi ); // insert into histogram q2 += step; } myHist->Output(filename); // write into file return myHist; } else { // read table and create histogram msg_Tracking()<<"HADRONS::Tau_Three_Pseudo::KS::CreateGHistogram : \n" <<" Read G(q2) from "<Extrapolate( s+p_hist->BinSize(), &val, 0 ); // shift +BinSize to get right value return val; } double ResonanceFlavour::IntegralG( double Q2, ResonanceFlavour res1, ResonanceFlavour res2, double beta, kf_code out ) { if (Q2==0.0) return 0.0; int Ns=500, Nt=500; // number of subintervals double sum (0.); double mpi2 = sqr( Flavour(out).HadMass() ); double s_max = Q2 + mpi2 - 2.*sqrt(Q2*mpi2); double s_min = 4.*mpi2; double ds = (s_max-s_min)/Ns; double t_max (0.); double t_min (0.); double dt (0.); double s (s_min), t, u; double V12, V22, V1V2; while ( s4.*ms && m_mass2>4.*ms) return( m_width*m_mass2/sqrt(s) * pow( (s-4.*ms)/(m_mass2-4.*ms), 0.5 ) ); return 0.; } double ResonanceFlavour::TwoBodyResonanceMassWidth( double s, double m ) { double ms=sqr(m); if (s>4.*ms && m_mass2>4.*ms) return( sqrt(s)*m_width*m_mass2/s * pow( (s-4.*ms)/(m_mass2-4.*ms), 1.5 ) ); return 0.; } double ResonanceFlavour::TwoBodyResonanceMassWidth( double s, double m1, double m2 ) { double threshold = sqr(m1+m2); double ms1 (m1*m1), ms2 (m2*m2); if (m_mass2>threshold && s>threshold) return( sqrt(s)*m_width*m_mass2/s * pow( m_mass2/s*Lambda(s,ms1,ms2)/Lambda(m_mass2,ms1,ms2), 1.5 ) ); return 0; } Complex ResonanceFlavour::BreitWigner( double s ) { double MG; if( m_running ) { MG = OffShellMassWidth(s); return Complex(m_mass2,0.)/Complex( m_mass2-s, -1.*MG ); } else { MG = m_mass * m_width; return Complex(m_mass2, -1.*MG)/Complex( m_mass2-s, -1.*MG ); } } Complex ResonanceFlavour::BreitWignerAlt( double s ) { double MG = m_mass * m_width; return Complex(1.,0.)/Complex( m_mass2-s, -1.*MG ); }