00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014 #include "EvtGenBase/EvtPatches.hh"
00015
00016
00017
00018
00019
00020
00021
00022 #include <stdlib.h>
00023 #include "EvtGenModels/EvtBtoXsgamma.hh"
00024 #include "EvtGenModels/EvtBtoXsgammaFlatEnergy.hh"
00025 #include "EvtGenBase/EvtRandom.hh"
00026 #include <fstream>
00027 using std::endl;
00028 using std::fstream;
00029
00030 EvtBtoXsgammaFlatEnergy::~EvtBtoXsgammaFlatEnergy(){
00031 }
00032
00033 void EvtBtoXsgammaFlatEnergy::init(int nArg, double* args){
00034
00035 if ((nArg) > 3 || (nArg > 1 && nArg <3)){
00036
00037 report(ERROR,"EvtGen") << "EvtBtoXsgamma generator model "
00038 << "EvtBtoXsgammaFlatEnergy expected "
00039 << "either 1(default config) or two arguments but found: "<<nArg<<endl;
00040 report(ERROR,"EvtGen") << "Will terminate execution!"<<endl;
00041 ::abort();
00042 }
00043 _mB0=5.2794;
00044 double mPi = 0.140;
00045 double mK = 0.494;
00046 if(nArg == 1){
00047 _eMin = 1.7;
00048
00049
00050 _eMax = (pow(_mB0,2)-pow(mPi+mK,2))/(2.0*_mB0);
00051 }else{
00052 _eMin=args[1];
00053 _eMax=args[2];
00054 }
00055 if (_eMax>(pow(_mB0,2)-pow(mPi+mK,2))/(2.0*_mB0)){
00056 report(ERROR,"EvtGen") << "Emax greater than Kinematic limit" << endl;
00057 report(ERROR,"EvtGen") << "Reset to the kinematic limit" << endl;
00058 report(ERROR,"EvtGen") << "(m_B**2-(m_pi+m_k)**2)/(2m_B)" << endl;
00059 _eMax = (pow(_mB0,2)-pow(mPi+mK,2))/(2.0*_mB0);
00060 }
00061 _eRange=_eMax-_eMin;
00062 }
00063
00064 double EvtBtoXsgammaFlatEnergy::GetMass( int Xscode ){
00065
00066 double eGamma = EvtRandom::Flat(_eRange)+_eMin;
00067 double mH = sqrt(pow(_mB0,2)-2.0*_mB0*eGamma);
00068 return mH;
00069 }
00070