00001
00002
00003
00004
00005 #include "stdlib.h"
00006 #include <iostream>
00007 #include "math.h"
00008 #include <map>
00009
00010 #include "StGenParticle.h"
00011
00012 void StGenParticle::Print(const char *opt) const
00013 {
00014 static int nCall=0; nCall++;
00015 std::cout << GetIdx() << " -" << std::flush;
00016 std::cout << " Ist=" << GetStatusCode() << std::flush;
00017 std::cout << " Pdg=" << GetPdgCode() << std::flush;
00018 std::cout << " Gea=" << GetGeaCode() << std::flush;
00019
00020 double V[3]; Vertex(V); std::cout << " Z=" << V[2] << std::flush;
00021
00022 if (GetNDaughters()) std::cout << "\tKids=" << GetNDaughters() << std::flush;
00023
00024 int moth1 = -1,moth2=-1;
00025 const StGenParticle *m = GetMother(0);
00026 if (m) moth1 = m->GetIdx();
00027 m = GetMother(1);
00028 if (m) moth2 = m->GetIdx();
00029
00030 if (moth1>=0 || moth2>=0) {
00031 std::cout << "\tMoth=(" << std::flush;
00032 if (moth1>=0) {std::cout << moth1 << std::flush;} else {std::cout << "_" << std::flush;}
00033 std::cout << " " << std::flush;
00034 if (moth2>=0) {std::cout << moth2 << std::flush;} else {std::cout << "_" << std::flush;}
00035 std::cout << ")" << std::flush;
00036 }
00037
00038 std::cout << std::endl;
00039 }
00040
00041 double StGenParticle::GetCalcMass() const
00042 {
00043 double p[4];
00044 Momentum(p);
00045 double m =p[3]*p[3]-(p[0]*p[0]+p[1]*p[1]+p[2]*p[2]);
00046 return (m>0) ? sqrt(m) : -sqrt(-m);
00047 }
00048
00049 double StGenParticle::R() const
00050 {
00051 double x[3]; Vertex(x); return sqrt(x[0]*x[0]+x[1]*x[1]);
00052
00053 }
00054
00055 double StGenParticle::Rho() const
00056 {
00057 double x[3]; Vertex(x); return sqrt(x[0]*x[0]+x[1]*x[1]+x[2]*x[2]);
00058
00059 }
00060
00061 double StGenParticle::P() const
00062 {
00063 double p[4]; Momentum(p); return sqrt(p[0]*p[0]+p[1]*p[1]+p[2]*p[2]);
00064
00065 }
00066
00067 double StGenParticle::Pt() const
00068 {
00069 double p[4]; Momentum(p); return sqrt(p[0]*p[0]+p[1]*p[1]);
00070
00071 }
00072
00073 double StGenParticle::Energy() const
00074 {
00075 double p[4]; Momentum(p); return p[3];
00076
00077 }
00078
00079 double StGenParticle::Eta() const
00080 {
00081 double p[4]; Momentum(p);
00082 double pmom = P();
00083 if (pmom > fabs(p[2])) return 0.5*log((pmom+p[2])/(pmom-p[2]));
00084 else return 1.e30;
00085 }
00086
00087 double StGenParticle::Phi() const
00088 {
00089 double p[4]; Momentum(p);
00090 return atan2(p[1],p[0]);
00091 }
00092
00093 double StGenParticle::Theta() const
00094 {
00095 double p[4]; Momentum(p);
00096 return acos(p[2]/P());
00097 }
00098 int StGenParticle::GetPdgCode() const { return StGenParticleMaster::Gea2Pdg(GetGeaCode());}
00099 int StGenParticle::GetGeaCode() const { return StGenParticleMaster::Pdg2Gea(GetPdgCode());}
00100
00101
00102
00103 void StGenParticleMaster::Print(const char *tit) const
00104 {
00105 if (!tit) tit = "";
00106 std::cout << "StGenParticleMaster::Print(" << tit << ")" << std::endl;
00107 const StGenParticle *p=0;
00108 for (int i=0;p=(*this)(i);i++) {
00109 p->Print();
00110 }
00111 std::cout << std::endl;
00112
00113 }
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124 enum {kNGEA = 50};
00125 static const int IDGEA[] = {
00126 0, 1, 8, 9, 7, 10, 11, 12, 16, 17,
00127 14, 13, 15, 25, 2, 3, 4, 4, 5, 6,
00128 18, 19, 20, 21, 22, 23, 24,
00129 26, 27, 28, 29, 30, 31, 32,
00130 33, 34, 35, 36, 37, 38, 39, 40, 41,
00131 42, 43, 44, 4, 4, 4, 4, 48, 48, 48,
00132 45, 46, 47, 49};
00133
00134 static const int IDPDG[] = {
00135 0, 22, 211, -211, 111, 130, 321, -321, 310, 221,
00136 2212, 2112,-2212,-2112, -11, 11, -12, 12, -13, 13,
00137 3122, 3222, 3212, 3112, 3322, 3312, 3334,
00138 -3122,-3222,-3212,-3112,-3322,-3312,-3334,
00139 -15, 15, 411, -411, 421, -421, 431, -431, 4122,
00140 24, -24, 23, -14, 14, -16, 16, 71, 72, 75,
00141 700201,700301,700202,700302};
00142
00143 int StGenParticleMaster::Gea2Pdg(int igea)
00144 {
00145
00146 static int PdgOfGea[kNGEA]={0};
00147
00148 if (! PdgOfGea[1]) {
00149 for (int i = 0;i<kNGEA; i++) { PdgOfGea[IDGEA[i]] = IDPDG[i]; }
00150 }
00151 if (igea >=kNGEA) return 0;
00152 return PdgOfGea[igea];
00153 }
00154
00155 double StGenParticleMaster::Gea2Mas(int igea)
00156 {
00157 static double mass[kNGEA]={
00158 0 ,0 ,0.00051 ,0.00051 ,0
00159 ,0.10566 ,0.10566 ,0.13498 ,0.13957 ,0.13957
00160 ,0.49767 ,0.4936 ,0.4936 ,0.93957 ,0.93827
00161 ,0.93827 ,0.49767 ,0.54745 ,1.11568 ,1.18937
00162 ,1.19255 ,1.19744 ,1.3149 ,1.3213 ,1.67245
00163 ,0.93957 ,1.11568 ,1.18937 ,1.19255 ,1.19744
00164 ,1.3149 ,1.3213 ,1.67245 ,1.777 ,1.777
00165 ,1.8693 ,1.8693 ,1.8645 ,1.8645 ,1.9685
00166 ,1.9685 ,2.2849 ,80.33 ,80.33 ,91.187
00167 ,0 ,0 ,0 ,0 ,0};
00168
00169 if (igea>=kNGEA) return 0.;
00170 return mass[igea];
00171 }
00172
00173
00174 int StGenParticleMaster::Pdg2Gea(int ipdg)
00175 {
00176 typedef std::map<int, int > pdgMap_t;
00177 static pdgMap_t pdgMap;
00178 static int once = 2009;
00179
00180 if (once) { once=0; for (int i=0;i<kNGEA;i++) {pdgMap[IDPDG[i]] = IDGEA[i];}}
00181
00182 pdgMap_t::iterator it = pdgMap.find(ipdg);
00183 if (it == pdgMap.end()) return 0;
00184 return (*it).second;
00185 }
