BOSS 8.0.0
BESIII Offline Software System
Loading...
Searching...
No Matches
EvtTrackGen Class Reference

#include <EvtTrackGen.hh>

Inheritance diagram for EvtTrackGen:

Public Member Functions

 EvtTrackGen ()
virtual ~EvtTrackGen ()
void getName (std::string &name)
EvtDecayBaseclone ()
void initProbMax ()
void init ()
void decay (EvtParticle *p)
Public Member Functions inherited from EvtDecayIncoherent
void makeDecay (EvtParticle *p)
virtual ~EvtDecayIncoherent ()
void setDaughterSpinDensity (int daughter)
int isDaughterSpinDensitySet (int daughter)
Public Member Functions inherited from EvtDecayBase
virtual std::string commandName ()
virtual void command (std::string cmd)
double getProbMax (double prob)
double resetProbMax (double prob)
 EvtDecayBase ()
virtual ~EvtDecayBase ()
virtual bool matchingDecay (const EvtDecayBase &other) const
EvtId getParentId ()
double getBranchingFraction ()
void disableCheckQ ()
void checkQ ()
int getNDaug ()
EvtIdgetDaugs ()
EvtId getDaug (int i)
int getNArg ()
int getPHOTOS ()
void setPHOTOS ()
void setVerbose ()
void setSummary ()
double * getArgs ()
std::string * getArgsStr ()
double getArg (int j)
std::string getArgStr (int j)
std::string getModelName ()
int getDSum ()
int summary ()
int verbose ()
void saveDecayInfo (EvtId ipar, int ndaug, EvtId *daug, int narg, std::vector< std::string > &args, std::string name, double brfr)
void printSummary ()
void setProbMax (double prbmx)
void noProbMax ()
void checkNArg (int a1, int a2=-1, int a3=-1, int a4=-1)
void checkNDaug (int d1, int d2=-1)
void checkSpinParent (EvtSpinType::spintype sp)
void checkSpinDaughter (int d1, EvtSpinType::spintype sp)
virtual int nRealDaughters ()

Additional Inherited Members

Static Public Member Functions inherited from EvtDecayBase
static void findMasses (EvtParticle *p, int ndaugs, EvtId daugs[10], double masses[10])
static void findMass (EvtParticle *p)
static double findMaxMass (EvtParticle *p)
Protected Member Functions inherited from EvtDecayBase
bool daugsDecayedByParentModel ()
Protected Attributes inherited from EvtDecayBase
bool _daugsDecayedByParentModel

Detailed Description

Definition at line 35 of file EvtTrackGen.hh.

Constructor & Destructor Documentation

◆ EvtTrackGen()

EvtTrackGen::EvtTrackGen ( )
inline

Definition at line 38 of file EvtTrackGen.hh.

38{}

Referenced by clone().

◆ ~EvtTrackGen()

EvtTrackGen::~EvtTrackGen ( )
virtual

Definition at line 25 of file EvtTrackGen.cc.

25{}

Member Function Documentation

◆ clone()

EvtDecayBase * EvtTrackGen::clone ( )
virtual

Implements EvtDecayBase.

Definition at line 29 of file EvtTrackGen.cc.

29{ return new EvtTrackGen; }

◆ decay()

void EvtTrackGen::decay ( EvtParticle * p)
virtual

Implements EvtDecayBase.

Definition at line 105 of file EvtTrackGen.cc.

105 {
106
107 double weight = p->initializePhaseSpace( getNDaug(), getDaugs() );
108 int rdm = (int)Evt.size() * EvtRandom::Flat( 0.0, 1.0 );
109 if ( Evt.size() == 0 )
110 {
111 std::cout << "EvtTrackGen: out of stored file record" << std::endl;
112 abort();
113 }
114 EvtVector4R ptot( 0, 0, 0, 0 );
115 for ( int i = 0; i < nParticles; i++ )
116 {
117 EvtParticle* daug = p->getDaug( i );
118 ptot += Evt[rdm][i];
119 daug->init( daug->getId(), Evt[rdm][i] );
120 }
121 // p->init(p->getId(),ptot);
122
123 // debugging
124 // std::cout<<p->getDaug(getNDaug()-1)->getP4()<<" =? "<<eParticle<<" "<<pxParticle<<"
125 // "<<pyParticle<<"
126 // "<<pzParticle<<std::endl;
127
128 return;
129}
*********Class see also m_nmax DOUBLE PRECISION m_MasPhot DOUBLE PRECISION m_phsu DOUBLE PRECISION m_Xenph DOUBLE PRECISION m_r2 DOUBLE PRECISION m_WtMass INTEGER m_nmax INTEGER m_Nevgen INTEGER m_IsFSR INTEGER m_MarTot *COMMON c_KarFin $ !Output file $ !Event serial number $ !alpha QED at Thomson limit $ !minimum energy at CMS for remooval $ !infrared dimensionless $ !dummy photon IR regulator $ !crude photon multiplicity enhancement factor *EVENT $ !MC crude volume of PhhSpace *Sfactors $ !YFS formfactor IR part only $ !YFS formfactor non IR finite part $ !mass weight
Definition KarFin.h:34
EvtId * getDaugs()
virtual void init(EvtId part_n, const EvtVector4R &p4)=0
EvtId getId() const
EvtParticle * getDaug(int i)
double initializePhaseSpace(int numdaughter, EvtId *daughters, double poleSize=-1., int whichTwo1=0, int whichTwo2=1)
static double Flat()
Definition EvtRandom.cc:69

◆ getName()

void EvtTrackGen::getName ( std::string & name)
virtual

Implements EvtDecayBase.

Definition at line 27 of file EvtTrackGen.cc.

27{ model_name = "TrackGen"; }

◆ init()

void EvtTrackGen::init ( )
virtual

Reimplemented from EvtDecayBase.

Definition at line 31 of file EvtTrackGen.cc.

31 {
32
33 // check that there are 0 arguments
34 checkNArg( 1 );
35 int idx = getArg( 0 );
36 if ( EvtGlobalSet::SV.size() == 0 )
37 {
38 std::cout << "No track data file is available! " << std::endl;
39 abort();
40 }
41 else { m_inputFileName = EvtGlobalSet::SV[idx]; }
42
43 m_inputFile.open( m_inputFileName.c_str() );
44 if ( !m_inputFile )
45 {
46 cout << "EvtTrackGen: PROBLEMS OPENING FILE " << m_inputFileName << endl;
47 exit( 0 );
48 }
49 // load the event P4
50 Evt.clear();
51 while ( !m_inputFile.eof() )
52 {
53 m_inputFile >> nParticles;
54 std::vector<EvtVector4R> vp4;
55 for ( int i = 0; i < nParticles; i++ )
56 {
57 EvtVector4R p4;
58 m_inputFile >> idParticles[i];
59 m_inputFile >> pxParticle;
60 p4.set( 1, pxParticle );
61 m_inputFile >> pyParticle;
62 p4.set( 2, pyParticle );
63 m_inputFile >> pzParticle;
64 p4.set( 3, pzParticle );
65 m_inputFile >> eParticle;
66 p4.set( 0, eParticle );
67 vp4.push_back( p4 );
68 }
69 Evt.push_back( vp4 );
70 }
71 Evt.pop_back();
72 // check daughters id
73 if ( nParticles != getNDaug() )
74 {
75 std::cout << "The number of daughters are not cosistent with that the data file"
76 << std::endl;
77 abort();
78 }
79 for ( int i = 0; i < nParticles; i++ )
80 {
81 EvtId pid = EvtPDL::evtIdFromStdHep( idParticles[i] );
82 if ( pid != getDaug( i ) )
83 {
84 std::cout << "The daughter particle pdg in your data file is not consistent with you "
85 "decay card."
86 << std::endl;
87 abort();
88 }
89 }
90 // debugging
91 /*
92 for(int i=0;i<Evt.size();i++){
93 std::cout<<"Event "<<i<<std::endl;
94 for(int j=0;j<nParticles;j++){
95 std::cout<<Evt[i][j].get(0)<<" "<<Evt[i][j].get(1)<<" "<<Evt[i][j].get(2)<<"
96 "<<Evt[i][j].get(3)<<std::endl;
97 }
98 }
99 */
100 //---------
101}
double getArg(int j)
void checkNArg(int a1, int a2=-1, int a3=-1, int a4=-1)
EvtId getDaug(int i)
static std::vector< std::string > SV
static EvtId evtIdFromStdHep(int stdhep)
Definition EvtPDL.cc:232
void set(int i, double d)

◆ initProbMax()

void EvtTrackGen::initProbMax ( )
virtual

Reimplemented from EvtDecayBase.

Definition at line 103 of file EvtTrackGen.cc.

103{ noProbMax(); }

The documentation for this class was generated from the following files: