00001 #include "Validation/EventGenerator/interface/TauDecay.h" 00002 #include "Validation/EventGenerator/interface/PdtPdgMini.h" 00003 00004 00005 #include <iomanip> 00006 #include <cstdlib> 00007 #include <iostream> 00008 00009 TauDecay::TauDecay(){ 00010 Reset(); 00011 } 00012 00013 TauDecay::~TauDecay(){ 00014 00015 } 00016 00017 void TauDecay::Reset(){ 00018 n_pi=0; 00019 n_pi0=0; 00020 n_K=0; 00021 n_K0L=0; 00022 n_K0S=0; 00023 n_gamma=0; 00024 n_nu=0; 00025 n_e=0; 00026 n_mu=0; 00027 n_a1=0; 00028 n_a10=0; 00029 n_rho=0; 00030 n_rho0=0; 00031 n_eta=0; 00032 n_omega=0; 00033 n_Kstar0=0; 00034 n_Kstar=0; 00035 unknown=0; 00036 } 00037 00038 bool TauDecay::isTauFinalStateParticle(int pdgid){ 00039 int id=abs(pdgid); 00040 if(id==PdtPdgMini::e_minus) return true; // e+- 00041 if(id==PdtPdgMini::nu_e) return true; // nu_e 00042 if(id==PdtPdgMini::mu_minus) return true; // mu+- 00043 if(id==PdtPdgMini::nu_mu) return true; // nu_mu 00044 if(id==PdtPdgMini::nu_tau) return true; // nu_tau 00045 if(id==PdtPdgMini::gamma) return true; // gamma happends in generator 00046 if(id==PdtPdgMini::pi0) return true; // pi0 00047 if(id==PdtPdgMini::pi_plus) return true; // pi+- 00048 if(id==PdtPdgMini::K_L0) return true; // K0L 00049 if(id==PdtPdgMini::K_S0) return true; // KS 00050 if(id==PdtPdgMini::K_plus) return true; // K+- 00051 return false; 00052 } 00053 00054 bool TauDecay::isTauParticleCounter(int pdgid){ 00055 int id=abs(pdgid); 00056 //count particles 00057 if(id==PdtPdgMini::pi_plus) { n_pi++; return true;} 00058 if(id==PdtPdgMini::pi0) { n_pi0++; return true;} 00059 if(id==PdtPdgMini::K_plus) { n_K++; return true;} 00060 if(id==PdtPdgMini::K_L0) { n_K0L++; return true;} 00061 if(id==PdtPdgMini::K_S0) { n_K0S++; return true;} 00062 if(id==PdtPdgMini::gamma) { n_gamma++; return true;} 00063 if(id==PdtPdgMini::nu_tau || 00064 id==PdtPdgMini::nu_e || 00065 id==PdtPdgMini::nu_mu) { n_nu++; return true;} 00066 if(id==PdtPdgMini::e_minus) { n_e++; return true;} 00067 if(id==PdtPdgMini::mu_minus) { n_mu++; return true;} 00068 return false; 00069 } 00070 00071 bool TauDecay::isTauResonanceCounter(int pdgid){ 00072 int id=abs(pdgid); 00073 //count resonances 00074 if(id==PdtPdgMini::a_1_plus) { n_a1++; return true;} 00075 if(id==PdtPdgMini::a_10) { n_a10++; return true;} 00076 if(id==PdtPdgMini::rho_plus) { n_rho++; return true;} 00077 if(id==PdtPdgMini::rho0) { n_rho0++; return true;} 00078 if(id==PdtPdgMini::eta) { n_eta++; return true;} 00079 if(id==PdtPdgMini::omega) { n_omega++; return true;} 00080 //if(id==PdtPdgMini::K_S0) { n_K0S++; return true;} 00081 if(id==PdtPdgMini::K_star0) { n_Kstar0++; return true;} 00082 if(id==PdtPdgMini::K_star_plus){ n_Kstar++; return true;} 00083 if(id==PdtPdgMini::W_plus) { return true;} 00084 unknown++; 00085 return false; 00086 } 00087 00088 void TauDecay::ClassifyDecayMode(unsigned int &JAK_ID,unsigned int &TauBitMask){ 00089 //Reset Bits 00090 JAK_ID=0; 00091 TauBitMask=0; 00092 // Classify according to JAK and TauDecayStructure 00094 // 00095 // Exlusive modes remove first 00096 // 00097 if(n_pi>=1 && n_pi0>=1 && n_nu==1 && n_eta==1){ // eta modes 00098 JAK_ID=JAK_ETAPIPI0; 00099 TauBitMask=OneProng; 00100 if(n_pi0==1)TauBitMask+=OnePi0; 00101 if(n_pi0==2)TauBitMask+=TwoPi0; 00102 if(n_pi0==3)TauBitMask+=ThreePi0; 00103 ClassifyDecayResonance(TauBitMask); 00104 return; 00105 } 00106 //JAK_K0BK0PI 00107 if ((n_K0S+n_K0L==2 && n_pi==1 && n_pi0==0 && n_K==0 && n_nu==1) || (n_pi>=1 && n_pi0==0 && n_K0L+n_K0S==2 && n_K==0 && n_nu==1) ){ 00108 JAK_ID=JAK_K0BK0PI; 00109 TauBitMask=OneProng; 00110 if(n_pi>1)TauBitMask+=KS0_to_pipi; 00111 ClassifyDecayResonance(TauBitMask); 00112 return; 00113 } 00114 //if(n_Kstar==1){ //JAK_KSTAR K0SPi && KPi0 00115 if(n_e==0 && n_mu==0 && n_pi>=1 && n_pi0==0 && n_K==0 && n_K0L+n_K0S==1 && n_nu==1){ 00116 JAK_ID=JAK_KSTAR; 00117 TauBitMask=OneProng; 00118 if(n_pi==3)TauBitMask+=KS0_to_pipi; 00119 ClassifyDecayResonance(TauBitMask); 00120 return; 00121 } 00122 if(n_e==0 && n_mu==0 && n_pi==0 && n_pi0==1 && n_K==1 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00123 JAK_ID=JAK_KSTAR; 00124 TauBitMask=OneProng; 00125 ClassifyDecayResonance(TauBitMask); 00126 return; 00127 } 00128 //} 00129 //JAK_PIK0PI0 00130 if(n_e==0 && n_mu==0 && n_pi>=1 && n_pi0==1 && n_K==0 && n_K0L+n_K0S==1 && n_nu==1){ 00131 JAK_ID=JAK_PIK0PI0; 00132 TauBitMask=OneProng; 00133 TauBitMask+=OnePi0; 00134 if(n_pi==3)TauBitMask+=KS0_to_pipi; 00135 ClassifyDecayResonance(TauBitMask); 00136 return; 00137 } 00138 //JAK_KK0B 00139 if(n_e==0 && n_mu==0 && n_pi0==0 && n_K==1 && n_K0L+n_K0S==1 && n_nu==1){ 00140 JAK_ID=JAK_KK0B; 00141 TauBitMask=OneProng; 00142 if(n_pi==2)TauBitMask+=KS0_to_pipi; 00143 return; 00144 } 00145 //JAK_ID=JAK_KK0BPI0 00146 if(n_e==0 && n_mu==0 && n_pi0==1 && n_K==1 && n_K0L+n_K0S==1 && n_nu==1){ 00147 JAK_ID=JAK_KK0BPI0; 00148 TauBitMask=OneProng; 00149 TauBitMask+=OnePi0; 00150 if(n_pi==2)TauBitMask+=KS0_to_pipi; 00151 ClassifyDecayResonance(TauBitMask); 00152 return; 00153 } 00154 00155 00156 //Safty handelling for exlusive modes 00157 if (n_K0L!=0){ 00158 std::cout << "Unknown mode with KL0: n_e " << n_e << " n_mu " << n_mu << " n_pi " << n_pi << " n_pi0 " << n_pi0 << " n_K " << n_K << " n_K0L " << n_K0L << " n_K0S " << n_K0S << " n_nu " << n_nu << " n_gamma " << n_gamma 00159 << std::endl; 00160 return; 00161 } 00162 if (n_K0S!=0){ 00163 std::cout << "Unknown mode with KS0: n_e " << n_e << " n_mu " << n_mu << " n_pi " << n_pi << " n_pi0 " << n_pi0 << " n_K " << n_K << " n_K0L " << n_K0L << " n_K0S " << n_K0S << " n_nu " << n_nu << " n_gamma " << n_gamma 00164 << std::endl; 00165 return; 00166 } 00167 if (n_eta!=0){ 00168 std::cout << "Unknown mode with eta: n_e " << n_e << " n_mu " << n_mu << " n_pi " << n_pi << " n_pi0 " << n_pi0 << " n_K " << n_K << " n_K0L " << n_K0L << " n_K0S " << n_K0S << " n_nu " << n_nu << " n_gamma " << n_gamma 00169 << std::endl; 00170 return; 00171 } 00172 00173 00174 00175 if(n_pi+n_K+n_e+n_mu==1)TauBitMask=OneProng; 00176 if(n_pi+n_K==3)TauBitMask=ThreeProng; 00177 if(n_pi+n_K==5)TauBitMask=FiveProng; 00178 if(n_pi0==1)TauBitMask+=OnePi0; 00179 if(n_pi0==2)TauBitMask+=TwoPi0; 00180 if(n_pi0==3)TauBitMask+=ThreePi0; 00181 ClassifyDecayResonance(TauBitMask); 00183 // 00184 // Standard modes 00185 // 00186 if(n_e==1 && n_mu==0 && n_pi==0 && n_pi0==0 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==2){ 00187 JAK_ID=JAK_ELECTRON; 00188 return; 00189 } 00190 if(n_e==0 && n_mu==1 && n_pi==0 && n_pi0==0 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==2){ 00191 JAK_ID=JAK_MUON; 00192 return; 00193 } 00194 if(n_e==0 && n_mu==0 && n_pi==1 && n_pi0==0 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00195 JAK_ID=JAK_PION; 00196 return; 00197 } 00198 if(n_e==0 && n_mu==0 && n_pi==1 && n_pi0==1 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==1 && n_rho==1){ 00199 JAK_ID=JAK_RHO_PIPI0; 00200 return; 00201 } 00202 if(n_e==0 && n_mu==0 && n_pi==1 && n_pi0==2 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00203 JAK_ID=JAK_A1_3PI; 00204 return; 00205 } 00206 if(n_e==0 && n_mu==0 && n_pi==3 && n_pi0==0 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00207 JAK_ID=JAK_A1_3PI; 00208 return; 00209 } 00210 if(n_e==0 && n_mu==0 && n_pi==0 && n_pi0==0 && n_K==1 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00211 JAK_ID=JAK_KAON; 00212 return; 00213 } 00214 if(n_e==0 && n_mu==0 && n_pi==3 && n_pi0==1 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00215 JAK_ID=JAK_3PIPI0; 00216 return; 00217 } 00218 if(n_e==0 && n_mu==0 && n_pi==1 && n_pi0==3 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00219 JAK_ID=JAK_PI3PI0; 00220 return; 00221 } 00222 if(n_e==0 && n_mu==0 && n_pi==3 && n_pi0==2 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00223 JAK_ID=JAK_3PI2PI0; 00224 return; 00225 } 00226 if(n_e==0 && n_mu==0 && n_pi==5 && n_pi0==0 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00227 JAK_ID=JAK_5PI; 00228 return; 00229 } 00230 if(n_e==0 && n_mu==0 && n_pi==5 && n_pi0==1 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00231 JAK_ID=JAK_5PIPI0; 00232 return; 00233 } 00234 if(n_e==0 && n_mu==0 && n_pi==3 && n_pi0==3 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00235 JAK_ID=JAK_3PI3PI0; 00236 return; 00237 } 00238 if(n_e==0 && n_mu==0 && n_pi==1 && n_pi0==0 && n_K==2 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00239 JAK_ID=JAK_KPIK; 00240 return; 00241 } 00242 if(n_e==0 && n_mu==0 && n_pi==0 && n_pi0==2 && n_K==1 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00243 JAK_ID=JAK_K2PI0; 00244 return; 00245 } 00246 if(n_e==0 && n_mu==0 && n_pi==2 && n_pi0==0 && n_K==1 && n_K0L==0 && n_K0S==0 && n_nu==1){ 00247 JAK_ID=JAK_KPIPI; 00248 return; 00249 } 00250 if(n_e==0 && n_mu==0 && n_pi==1 && n_pi0==1 && n_K==0 && n_K0L==0 && n_K0S==0 && n_nu==1 && n_gamma>=1 && n_rho==0){ 00251 JAK_ID=JAK_PIPI0GAM; 00252 return; 00253 } 00254 std::cout << "Tau Mode not found: n_e " << n_e << " n_mu " << n_mu << " n_pi " << n_pi << " n_pi0 " << n_pi0 << " n_K " << n_K << " n_K0L " << n_K0L << " n_K0S " << n_K0S << " n_nu " << n_nu << " n_gamma " << n_gamma << std::endl; 00255 JAK_ID=JAK_UNKNOWN; 00256 } 00257 00258 void TauDecay::ClassifyDecayResonance(unsigned int &TauBitMask){ 00259 // Add Resonance info to TauBitMask 00260 if(n_a1>0) TauBitMask+=Res_a1_pm; 00261 if(n_a10>0) TauBitMask+=Res_a1_0; 00262 if(n_rho>0) TauBitMask+=Res_rho_pm; 00263 if(n_rho0>0) TauBitMask+=Res_rho_0; 00264 if(n_eta>0) TauBitMask+=Res_eta; 00265 if(n_omega>0) TauBitMask+=Res_omega; 00266 if(n_Kstar>0) TauBitMask+=Res_Kstar_pm; 00267 if(n_Kstar0>0) TauBitMask+=Res_Kstar_0; 00268 }