ATLAS_2012_I943401.cc
Go to the documentation of this file.
00001 // -*- C++ -*- 00002 #include "Rivet/Analysis.hh" 00003 #include "Rivet/Tools/BinnedHistogram.hh" 00004 #include "Rivet/Projections/FinalState.hh" 00005 #include "Rivet/Projections/ChargedFinalState.hh" 00006 #include "Rivet/Projections/VisibleFinalState.hh" 00007 #include "Rivet/Projections/IdentifiedFinalState.hh" 00008 #include "Rivet/Projections/FastJets.hh" 00009 #include "Rivet/Projections/VetoedFinalState.hh" 00010 00011 namespace Rivet { 00012 00013 00014 class ATLAS_2012_I943401 : public Analysis { 00015 public: 00016 00017 /// @name Constructors etc. 00018 //@{ 00019 00020 /// Constructor 00021 00022 ATLAS_2012_I943401() 00023 : Analysis("ATLAS_2012_I943401") 00024 { } 00025 00026 //@} 00027 00028 00029 public: 00030 00031 /// @name Analysis methods 00032 //@{ 00033 00034 /// Book histograms and initialise projections before the run 00035 void init() { 00036 00037 // projection to find the electrons 00038 std::vector<std::pair<double, double> > eta_e; 00039 eta_e.push_back(make_pair(-2.47,2.47)); 00040 IdentifiedFinalState elecs(eta_e, 20.0*GeV); 00041 elecs.acceptIdPair(PID::ELECTRON); 00042 addProjection(elecs, "elecs"); 00043 00044 // projection to find the muons 00045 std::vector<std::pair<double, double> > eta_m; 00046 eta_m.push_back(make_pair(-2.4,2.4)); 00047 IdentifiedFinalState muons(eta_m, 10.0*GeV); 00048 muons.acceptIdPair(PID::MUON); 00049 addProjection(muons, "muons"); 00050 00051 // jet finder 00052 VetoedFinalState vfs; 00053 vfs.addVetoPairId(PID::MUON); 00054 addProjection(FastJets(vfs, FastJets::ANTIKT, 0.4), 00055 "AntiKtJets04"); 00056 00057 // all tracks (to do deltaR with leptons) 00058 addProjection(ChargedFinalState(-3.0,3.0,0.5*GeV),"cfs"); 00059 00060 // for pTmiss 00061 addProjection(VisibleFinalState(-4.5,4.5),"vfs"); 00062 00063 // book histograms 00064 00065 // counts in signal regions 00066 _count_OS_SR1 = bookHisto1D("count_OS_SR1", 1, 0., 1.); 00067 _count_OS_SR2 = bookHisto1D("count_OS_SR2", 1, 0., 1.); 00068 _count_OS_SR3 = bookHisto1D("count_OS_SR3", 1, 0., 1.); 00069 _count_SS_SR1 = bookHisto1D("count_SS_SR1", 1, 0., 1.); 00070 _count_SS_SR2 = bookHisto1D("count_SS_SR2", 1, 0., 1.); 00071 _count_FS_SR1 = bookHisto1D("count_FS_SR1", 1, 0., 1.); 00072 _count_FS_SR2 = bookHisto1D("count_FS_SR2", 1, 0., 1.); 00073 _count_FS_SR3 = bookHisto1D("count_FS_SR3", 1, 0., 1.); 00074 00075 // histograms from paper 00076 00077 _hist_mll_SS_D = bookHisto1D( 1,1,1); 00078 _hist_mll_SS_B = bookHisto1D( 1,1,2); 00079 _hist_eTmiss_SS_D = bookHisto1D( 2,1,1); 00080 _hist_eTmiss_SS_B = bookHisto1D( 2,1,2); 00081 _hist_mll_SS_2Jet_D = bookHisto1D( 3,1,1); 00082 _hist_mll_SS_2Jet_B = bookHisto1D( 3,1,2); 00083 _hist_njet_SS_D = bookHisto1D( 5,1,1); 00084 _hist_njet_SS_B = bookHisto1D( 5,1,2); 00085 _hist_pT_j1_SS_D = bookHisto1D( 6,1,1); 00086 _hist_pT_j1_SS_B = bookHisto1D( 6,1,2); 00087 _hist_pT_j2_SS_D = bookHisto1D( 7,1,1); 00088 _hist_pT_j2_SS_B = bookHisto1D( 7,1,2); 00089 _hist_pT_l1_SS_D = bookHisto1D( 8,1,1); 00090 _hist_pT_l1_SS_B = bookHisto1D( 8,1,2); 00091 _hist_pT_l2_SS_D = bookHisto1D( 9,1,1); 00092 _hist_pT_l2_SS_B = bookHisto1D( 9,1,2); 00093 _hist_mll_OS_D = bookHisto1D(10,1,1); 00094 _hist_mll_OS_B = bookHisto1D(10,1,2); 00095 _hist_eTmiss_OS_D = bookHisto1D(11,1,1); 00096 _hist_eTmiss_OS_B = bookHisto1D(11,1,2); 00097 _hist_eTmiss_3Jet_OS_D = bookHisto1D(12,1,1); 00098 _hist_eTmiss_3Jet_OS_B = bookHisto1D(12,1,2); 00099 _hist_eTmiss_4Jet_OS_D = bookHisto1D(13,1,1); 00100 _hist_eTmiss_4Jet_OS_B = bookHisto1D(13,1,2); 00101 _hist_njet_OS_D = bookHisto1D(14,1,1); 00102 _hist_njet_OS_B = bookHisto1D(14,1,2); 00103 _hist_pT_j1_OS_D = bookHisto1D(15,1,1); 00104 _hist_pT_j1_OS_B = bookHisto1D(15,1,2); 00105 _hist_pT_j2_OS_D = bookHisto1D(16,1,1); 00106 _hist_pT_j2_OS_B = bookHisto1D(16,1,2); 00107 _hist_pT_l1_OS_D = bookHisto1D(17,1,1); 00108 _hist_pT_l1_OS_B = bookHisto1D(17,1,2); 00109 _hist_pT_l2_OS_D = bookHisto1D(18,1,1); 00110 _hist_pT_l2_OS_B = bookHisto1D(18,1,2); 00111 //???? 00112 // <dataPointSet name="d04-x01-y01" dimension="2" path="/REF/ATLAS_2011_I943401" title="EVENTS/10 GEV" > 00113 // <dataPointSet name="d04-x01-y02" dimension="2" path="/REF/ATLAS_2011_I943401" title="EVENTS/10 GEV" > 00114 } 00115 00116 /// Perform the event analysis 00117 void analyze(const Event& event) { 00118 // event weight 00119 const double weight = event.weight(); 00120 00121 // get the jet candidates 00122 Jets cand_jets; 00123 foreach (const Jet& jet, 00124 applyProjection<FastJets>(event, "AntiKtJets04").jetsByPt(20.0*GeV) ) { 00125 if ( fabs( jet.eta() ) < 2.8 ) { 00126 cand_jets.push_back(jet); 00127 } 00128 } 00129 00130 // electron candidates 00131 Particles cand_e = 00132 applyProjection<IdentifiedFinalState>(event, "elecs").particlesByPt(); 00133 00134 // Discard jets that overlap with electrons 00135 Jets recon_jets; 00136 foreach ( const Jet& jet, cand_jets ) { 00137 bool away_from_e = true; 00138 foreach ( const Particle & e, cand_e ) { 00139 if ( deltaR(e.momentum(),jet.momentum()) <= 0.2 ) { 00140 away_from_e = false; 00141 break; 00142 } 00143 } 00144 if ( away_from_e ) recon_jets.push_back( jet ); 00145 } 00146 // get the charged tracks for isolation 00147 Particles chg_tracks = 00148 applyProjection<ChargedFinalState>(event, "cfs").particles(); 00149 00150 // Reconstructed electrons 00151 Particles recon_e; 00152 foreach ( const Particle & e, cand_e ) { 00153 // check not near a jet 00154 bool e_near_jet = false; 00155 foreach ( const Jet& jet, recon_jets ) { 00156 if ( deltaR(e.momentum(),jet.momentum()) < 0.4 ) { 00157 e_near_jet = true; 00158 break; 00159 } 00160 } 00161 if ( e_near_jet ) continue; 00162 // check the isolation 00163 double pTinCone = -e.pT(); 00164 foreach ( const Particle & track, chg_tracks ) { 00165 if ( deltaR(e.momentum(),track.momentum()) < 0.2 ) 00166 pTinCone += track.pT(); 00167 } 00168 if ( pTinCone < 0.1*e.momentum().perp() ) 00169 recon_e.push_back(e); 00170 } 00171 00172 // Reconstructed Muons 00173 Particles recon_mu; 00174 Particles cand_mu = 00175 applyProjection<IdentifiedFinalState>(event,"muons").particlesByPt(); 00176 foreach ( const Particle & mu, cand_mu ) { 00177 // check not near a jet 00178 bool mu_near_jet = false; 00179 foreach ( const Jet& jet, recon_jets ) { 00180 if ( deltaR(mu.momentum(),jet.momentum()) < 0.4 ) { 00181 mu_near_jet = true; 00182 break; 00183 } 00184 } 00185 if ( mu_near_jet ) continue; 00186 // isolation 00187 double pTinCone = -mu.pT(); 00188 foreach ( const Particle & track, chg_tracks ) { 00189 if ( deltaR(mu.momentum(),track.momentum()) < 0.2 ) 00190 pTinCone += track.pT(); 00191 } 00192 if ( pTinCone < 1.8*GeV ) 00193 recon_mu.push_back(mu); 00194 } 00195 00196 // pTmiss 00197 Particles vfs_particles 00198 = applyProjection<VisibleFinalState>(event, "vfs").particles(); 00199 FourMomentum pTmiss; 00200 foreach ( const Particle & p, vfs_particles ) { 00201 pTmiss -= p.momentum(); 00202 } 00203 double eTmiss = pTmiss.pT(); 00204 00205 // ATLAS calo problem 00206 if(rand()/static_cast<double>(RAND_MAX)<=0.42) { 00207 foreach ( const Particle & e, recon_e ) { 00208 double eta = e.eta(); 00209 double phi = e.momentum().azimuthalAngle(MINUSPI_PLUSPI); 00210 if(eta>-0.1&&eta<1.5&&phi>-0.9&&phi<-0.5) 00211 vetoEvent; 00212 } 00213 foreach ( const Jet & jet, recon_jets ) { 00214 double eta = jet.rapidity(); 00215 double phi = jet.momentum().azimuthalAngle(MINUSPI_PLUSPI); 00216 if(jet.momentum().perp()>40 && eta>-0.1&&eta<1.5&&phi>-0.9&&phi<-0.5) 00217 vetoEvent; 00218 } 00219 } 00220 00221 // Exactly two leptons for each event 00222 if ( recon_mu.size() + recon_e.size() != 2) 00223 vetoEvent; 00224 // two electrons highest pT > 25 00225 Particles recon_leptons; 00226 if(recon_e.size()==2&&recon_e[0].momentum().perp()>25.) { 00227 recon_leptons = recon_e; 00228 } 00229 // two muons highest pT > 20 00230 else if(recon_mu.size()==2&&recon_mu[0].momentum().perp()>20.) { 00231 recon_leptons = recon_mu; 00232 } 00233 else if(recon_e.size()==1 && recon_mu.size()==1 && 00234 (recon_e[0].momentum().perp()>25. ||recon_mu[0].momentum().perp()>20. )) { 00235 if(recon_mu[0].momentum().perp()<recon_e[0].momentum().perp()) { 00236 recon_leptons.push_back(recon_e [0]); 00237 recon_leptons.push_back(recon_mu[0]); 00238 } 00239 else { 00240 recon_leptons.push_back(recon_mu[0]); 00241 recon_leptons.push_back(recon_e [0]); 00242 } 00243 } 00244 // fails trigger 00245 else 00246 vetoEvent; 00247 00248 double mll = (recon_leptons[0].momentum()+recon_leptons[1].momentum()).mass(); 00249 // lepton pair mass > 12. 00250 if(mll < 12.) vetoEvent; 00251 00252 // same sign or opposite sign event 00253 int sign = recon_leptons[0].pdgId()*recon_leptons[1].pdgId(); 00254 00255 // same sign leptons 00256 if(sign>0) { 00257 _hist_mll_SS_D ->fill(mll ,weight); 00258 _hist_mll_SS_B ->fill(mll ,weight); 00259 _hist_eTmiss_SS_D->fill(eTmiss,weight); 00260 _hist_eTmiss_SS_B->fill(eTmiss,weight); 00261 if(recon_jets.size()>=2) { 00262 _hist_mll_SS_2Jet_D ->fill(mll ,weight); 00263 _hist_mll_SS_2Jet_B ->fill(mll ,weight); 00264 } 00265 _hist_njet_SS_D ->fill(recon_jets.size(),weight); 00266 _hist_njet_SS_B ->fill(recon_jets.size(),weight); 00267 if(!recon_jets.empty()) { 00268 _hist_pT_j1_SS_D->fill(recon_jets[0].momentum().perp(),weight); 00269 _hist_pT_j1_SS_B->fill(recon_jets[0].momentum().perp(),weight); 00270 } 00271 if(recon_jets.size()>2) { 00272 _hist_pT_j2_SS_D->fill(recon_jets[1].momentum().perp(),weight); 00273 _hist_pT_j2_SS_B->fill(recon_jets[1].momentum().perp(),weight); 00274 } 00275 _hist_pT_l1_SS_D->fill(recon_leptons[0].momentum().perp(),weight); 00276 _hist_pT_l1_SS_B->fill(recon_leptons[0].momentum().perp(),weight); 00277 _hist_pT_l2_SS_D->fill(recon_leptons[1].momentum().perp(),weight); 00278 _hist_pT_l2_SS_B->fill(recon_leptons[1].momentum().perp(),weight); 00279 // SS-SR1 00280 if(eTmiss>100.) { 00281 _count_SS_SR1->fill(0.5,weight); 00282 } 00283 // SS-SR2 00284 if(eTmiss>80. && recon_jets.size()>=2 && 00285 recon_jets[1].momentum().perp()>50.) { 00286 _count_SS_SR2->fill(0.5,weight); 00287 } 00288 } 00289 // opposite sign 00290 else { 00291 _hist_mll_OS_D->fill(mll ,weight); 00292 _hist_mll_OS_B->fill(mll ,weight); 00293 _hist_eTmiss_OS_D->fill(eTmiss,weight); 00294 _hist_eTmiss_OS_B->fill(eTmiss,weight); 00295 if(recon_jets.size()>=3){ 00296 _hist_eTmiss_3Jet_OS_D->fill(eTmiss,weight); 00297 _hist_eTmiss_3Jet_OS_B->fill(eTmiss,weight); 00298 } 00299 if(recon_jets.size()>=4){ 00300 _hist_eTmiss_4Jet_OS_D->fill(eTmiss,weight); 00301 _hist_eTmiss_4Jet_OS_B->fill(eTmiss,weight); 00302 } 00303 _hist_njet_OS_D->fill(recon_jets.size(),weight); 00304 _hist_njet_OS_B->fill(recon_jets.size(),weight); 00305 if(!recon_jets.empty()) { 00306 _hist_pT_j1_OS_D->fill(recon_jets[0].momentum().perp(),weight); 00307 _hist_pT_j1_OS_B->fill(recon_jets[0].momentum().perp(),weight); 00308 } 00309 if(recon_jets.size()>2) { 00310 _hist_pT_j2_OS_D->fill(recon_jets[1].momentum().perp(),weight); 00311 _hist_pT_j2_OS_B->fill(recon_jets[1].momentum().perp(),weight); 00312 } 00313 _hist_pT_l1_OS_D->fill(recon_leptons[0].momentum().perp(),weight); 00314 _hist_pT_l1_OS_B->fill(recon_leptons[0].momentum().perp(),weight); 00315 _hist_pT_l2_OS_D->fill(recon_leptons[1].momentum().perp(),weight); 00316 _hist_pT_l2_OS_B->fill(recon_leptons[1].momentum().perp(),weight); 00317 // different signal regions 00318 // OS-SR1 00319 if(eTmiss>250.) { 00320 _count_OS_SR1->fill(0.5,weight); 00321 } 00322 // OS-SR2 00323 if(eTmiss>220. && recon_jets.size()>=3 && 00324 recon_jets[0].momentum().perp()>80. && 00325 recon_jets[2].momentum().perp()>40.) { 00326 _count_OS_SR2->fill(0.5,weight); 00327 } 00328 // OS-SR3 00329 if(eTmiss>100. && recon_jets.size()>=4 && 00330 recon_jets[0].momentum().perp()>100. && 00331 recon_jets[3].momentum().perp()>70.) { 00332 _count_OS_SR3->fill(0.5,weight); 00333 } 00334 // same flavour analysis 00335 static const double beta = 0.75; 00336 static const double tau_e = 0.96; 00337 static const double tau_mu = 0.816; 00338 double fs_weight = weight; 00339 if(abs(recon_leptons[0].pdgId())==PID::ELECTRON && abs(recon_leptons[1].pdgId())==PID::ELECTRON) { 00340 fs_weight /= beta*(1.-sqr(1.-tau_e)); 00341 } 00342 else if(abs(recon_leptons[0].pdgId())==PID::MUON && abs(recon_leptons[1].pdgId())==PID::MUON) { 00343 fs_weight *= beta/(1.-sqr(1.-tau_mu)); 00344 } 00345 else { 00346 fs_weight /= -(1.-(1.-tau_e)*(1.-tau_mu)); 00347 } 00348 // FS-SR1 00349 if(eTmiss>80.&& (mll<80.||mll>100.)) { 00350 _count_FS_SR1->fill(0.5,fs_weight); 00351 } 00352 // FS-SR2 00353 if(eTmiss>80.&&recon_jets.size()>=2) { 00354 _count_FS_SR2->fill(0.5,fs_weight); 00355 } 00356 // FS-SR3 00357 if(eTmiss>250.) { 00358 _count_FS_SR3->fill(0.5,fs_weight); 00359 } 00360 } 00361 } 00362 00363 //@} 00364 00365 00366 void finalize() { 00367 00368 double norm = crossSection()/femtobarn*1.04/sumOfWeights(); 00369 // event counts 00370 scale(_count_OS_SR1,norm); 00371 scale(_count_OS_SR2,norm); 00372 scale(_count_OS_SR3,norm); 00373 scale(_count_SS_SR1,norm); 00374 scale(_count_SS_SR2,norm); 00375 scale(_count_FS_SR1,norm); 00376 scale(_count_FS_SR2,norm); 00377 scale(_count_FS_SR3,norm); 00378 // histograms 00379 scale(_hist_mll_SS_D ,norm*20.); 00380 scale(_hist_mll_SS_B ,norm*20.); 00381 scale(_hist_eTmiss_SS_D ,norm*20.); 00382 scale(_hist_eTmiss_SS_B ,norm*20.); 00383 scale(_hist_mll_SS_2Jet_D,norm*50.); 00384 scale(_hist_mll_SS_2Jet_B,norm*50.); 00385 scale(_hist_njet_SS_D ,norm ); 00386 scale(_hist_njet_SS_B ,norm ); 00387 scale(_hist_pT_j1_SS_D ,norm*20.); 00388 scale(_hist_pT_j1_SS_B ,norm*20.); 00389 scale(_hist_pT_j2_SS_D ,norm*20.); 00390 scale(_hist_pT_j2_SS_B ,norm*20.); 00391 scale(_hist_pT_l1_SS_D ,norm*5. ); 00392 scale(_hist_pT_l1_SS_B ,norm*5. ); 00393 scale(_hist_pT_l2_SS_D ,norm*5. ); 00394 scale(_hist_pT_l2_SS_B ,norm*5. ); 00395 00396 scale(_hist_mll_OS_D ,norm*10.); 00397 scale(_hist_mll_OS_B ,norm*10.); 00398 scale(_hist_eTmiss_OS_D ,norm*10.); 00399 scale(_hist_eTmiss_OS_B ,norm*10.); 00400 scale(_hist_eTmiss_3Jet_OS_D,norm*10.); 00401 scale(_hist_eTmiss_3Jet_OS_B,norm*10.); 00402 scale(_hist_eTmiss_4Jet_OS_D,norm*10.); 00403 scale(_hist_eTmiss_4Jet_OS_B,norm*10.); 00404 scale(_hist_njet_OS_D ,norm ); 00405 scale(_hist_njet_OS_B ,norm ); 00406 scale(_hist_pT_j1_OS_D ,norm*20.); 00407 scale(_hist_pT_j1_OS_B ,norm*20.); 00408 scale(_hist_pT_j2_OS_D ,norm*20.); 00409 scale(_hist_pT_j2_OS_B ,norm*20.); 00410 scale(_hist_pT_l1_OS_D ,norm*20.); 00411 scale(_hist_pT_l1_OS_B ,norm*20.); 00412 scale(_hist_pT_l2_OS_D ,norm*20.); 00413 scale(_hist_pT_l2_OS_B ,norm*20.); 00414 } 00415 00416 private: 00417 00418 /// @name Histograms 00419 //@{ 00420 Histo1DPtr _count_OS_SR1; 00421 Histo1DPtr _count_OS_SR2; 00422 Histo1DPtr _count_OS_SR3; 00423 Histo1DPtr _count_SS_SR1; 00424 Histo1DPtr _count_SS_SR2; 00425 Histo1DPtr _count_FS_SR1; 00426 Histo1DPtr _count_FS_SR2; 00427 Histo1DPtr _count_FS_SR3; 00428 00429 Histo1DPtr _hist_mll_SS_D; 00430 Histo1DPtr _hist_mll_SS_B; 00431 Histo1DPtr _hist_eTmiss_SS_D; 00432 Histo1DPtr _hist_eTmiss_SS_B; 00433 Histo1DPtr _hist_mll_SS_2Jet_D; 00434 Histo1DPtr _hist_mll_SS_2Jet_B; 00435 Histo1DPtr _hist_njet_SS_D; 00436 Histo1DPtr _hist_njet_SS_B; 00437 Histo1DPtr _hist_pT_j1_SS_D; 00438 Histo1DPtr _hist_pT_j1_SS_B; 00439 Histo1DPtr _hist_pT_j2_SS_D; 00440 Histo1DPtr _hist_pT_j2_SS_B; 00441 Histo1DPtr _hist_pT_l1_SS_D; 00442 Histo1DPtr _hist_pT_l1_SS_B; 00443 Histo1DPtr _hist_pT_l2_SS_D; 00444 Histo1DPtr _hist_pT_l2_SS_B; 00445 00446 Histo1DPtr _hist_mll_OS_D; 00447 Histo1DPtr _hist_mll_OS_B; 00448 Histo1DPtr _hist_eTmiss_OS_D; 00449 Histo1DPtr _hist_eTmiss_OS_B; 00450 Histo1DPtr _hist_eTmiss_3Jet_OS_D; 00451 Histo1DPtr _hist_eTmiss_3Jet_OS_B; 00452 Histo1DPtr _hist_eTmiss_4Jet_OS_D; 00453 Histo1DPtr _hist_eTmiss_4Jet_OS_B; 00454 Histo1DPtr _hist_njet_OS_D ; 00455 Histo1DPtr _hist_njet_OS_B ; 00456 Histo1DPtr _hist_pT_j1_OS_D; 00457 Histo1DPtr _hist_pT_j1_OS_B; 00458 Histo1DPtr _hist_pT_j2_OS_D; 00459 Histo1DPtr _hist_pT_j2_OS_B; 00460 Histo1DPtr _hist_pT_l1_OS_D; 00461 Histo1DPtr _hist_pT_l1_OS_B; 00462 Histo1DPtr _hist_pT_l2_OS_D; 00463 Histo1DPtr _hist_pT_l2_OS_B; 00464 //@} 00465 }; 00466 00467 // The hook for the plugin system 00468 DECLARE_RIVET_PLUGIN(ATLAS_2012_I943401); 00469 00470 } Generated on Fri Oct 25 2013 12:41:44 for The Rivet MC analysis system by 
1.7.6.1
|
