00001
00002 #include "Rivet/Analysis.hh"
00003 #include "Rivet/RivetAIDA.hh"
00004 #include "Rivet/Tools/Logging.hh"
00005 #include "Rivet/Projections/FastJets.hh"
00006 #include "Rivet/Projections/TotalVisibleMomentum.hh"
00007
00008 namespace Rivet {
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025 class D0_2004_S5992206 : public Analysis {
00026
00027 public:
00028
00029
00030
00031
00032
00033 D0_2004_S5992206() : Analysis("D0_2004_S5992206")
00034 {
00035 setBeams(PROTON, ANTIPROTON);
00036 }
00037
00038
00039
00040
00041
00042
00043
00044 void init() {
00045
00046 const FinalState fs(-3.0, 3.0);
00047 addProjection(fs, "FS");
00048 addProjection(FastJets(FinalState(), FastJets::D0ILCONE, 0.7), "Jets");
00049 addProjection(TotalVisibleMomentum(fs), "CalMET");
00050
00051
00052 VetoedFinalState vfs(fs);
00053 vfs.vetoNeutrinos();
00054 vfs.addVetoPairDetail(MUON, 1.0, MAXDOUBLE);
00055 addProjection(vfs, "VFS");
00056
00057
00058 _histJetAzimuth_pTmax75_100 = bookHistogram1D(1, 2, 1);
00059 _histJetAzimuth_pTmax100_130 = bookHistogram1D(2, 2, 1);
00060 _histJetAzimuth_pTmax130_180 = bookHistogram1D(3, 2, 1);
00061 _histJetAzimuth_pTmax180_ = bookHistogram1D(4, 2, 1);
00062 }
00063
00064
00065
00066 void analyze(const Event & event) {
00067
00068
00069 const JetAlg& jetpro = applyProjection<JetAlg>(event, "Jets");
00070 getLog() << Log::DEBUG << "Jet multiplicity before any pT cut = " << jetpro.size() << endl;
00071
00072 const Jets jets = jetpro.jetsByPt(40.0*GeV);
00073 if (jets.size() >= 2) {
00074 getLog() << Log::DEBUG << "Jet multiplicity after pT > 40 GeV cut = " << jets.size() << endl;
00075 } else {
00076 vetoEvent;
00077 }
00078 const double rap1 = jets[0].momentum().rapidity();
00079 const double rap2 = jets[1].momentum().rapidity();
00080 if (fabs(rap1) > 0.5 || fabs(rap2) > 0.5) {
00081 vetoEvent;
00082 }
00083 getLog() << Log::DEBUG << "Jet eta and pT requirements fulfilled" << endl;
00084 const double pT1 = jets[0].momentum().pT();
00085
00086 const TotalVisibleMomentum& caloMissEt = applyProjection<TotalVisibleMomentum>(event, "CalMET");
00087 getLog() << Log::DEBUG << "Missing Et = " << caloMissEt.momentum().pT()/GeV << endl;
00088 if (caloMissEt.momentum().pT() > 0.7*pT1) {
00089 vetoEvent;
00090 }
00091
00092 if (pT1/GeV >= 75.0) {
00093 const double weight = event.weight();
00094 const double dphi = deltaPhi(jets[0].momentum().phi(), jets[1].momentum().phi());
00095 if (inRange(pT1/GeV, 75.0, 100.0)) {
00096 _histJetAzimuth_pTmax75_100->fill(dphi, weight);
00097 } else if (inRange(pT1/GeV, 100.0, 130.0)) {
00098 _histJetAzimuth_pTmax100_130->fill(dphi, weight);
00099 } else if (inRange(pT1/GeV, 130.0, 180.0)) {
00100 _histJetAzimuth_pTmax130_180->fill(dphi, weight);
00101 } else if (pT1/GeV > 180.0) {
00102 _histJetAzimuth_pTmax180_->fill(dphi, weight);
00103 }
00104 }
00105
00106 }
00107
00108
00109
00110 void finalize() {
00111
00112 normalize(_histJetAzimuth_pTmax75_100);
00113 normalize(_histJetAzimuth_pTmax100_130);
00114 normalize(_histJetAzimuth_pTmax130_180);
00115 normalize(_histJetAzimuth_pTmax180_);
00116 }
00117
00118
00119
00120
00121 private:
00122
00123
00124
00125 AIDA::IHistogram1D* _histJetAzimuth_pTmax75_100;
00126 AIDA::IHistogram1D* _histJetAzimuth_pTmax100_130;
00127 AIDA::IHistogram1D* _histJetAzimuth_pTmax130_180;
00128 AIDA::IHistogram1D* _histJetAzimuth_pTmax180_;
00129
00130
00131 };
00132
00133
00134
00135
00136 AnalysisBuilder<D0_2004_S5992206> plugin_D0_2004_S5992206;
00137
00138 }