MC_ZZJETS.cc

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 #include "Rivet/Analyses/MC_JetAnalysis.hh"
00003 #include "Rivet/Projections/ZFinder.hh"
00004 #include "Rivet/Projections/FastJets.hh"
00005 #include "Rivet/Projections/VetoedFinalState.hh"
00006 #include "Rivet/Tools/Logging.hh"
00007 #include "Rivet/Tools/ParticleIdUtils.hh"
00008 #include "Rivet/RivetAIDA.hh"
00009 
00010 namespace Rivet {
00011 
00012   /// @brief MC validation analysis for Z[ee]Z[mumu] + jets events
00013   class MC_ZZJETS : public MC_JetAnalysis {
00014   public:
00015 
00016     /// Default constructor
00017     MC_ZZJETS()
00018       : MC_JetAnalysis("MC_ZZJETS", 4, "Jets")
00019     {
00020       setNeedsCrossSection(true);
00021     }
00022 
00023 
00024     /// @name Analysis methods
00025     //@{
00026 
00027     /// Book histograms
00028     void init() {
00029       ZFinder zeefinder(-3.5, 3.5, 25.0*GeV, ELECTRON, 65.0*GeV, 115.0*GeV, 0.2);
00030       addProjection(zeefinder, "ZeeFinder");
00031       ZFinder zmmfinder(-3.5, 3.5, 25.0*GeV, MUON, 65.0*GeV, 115.0*GeV, 0.2);
00032       addProjection(zmmfinder, "ZmmFinder");
00033       VetoedFinalState jetinput;
00034       jetinput
00035           .addVetoOnThisFinalState(zeefinder.constituentsFinalState())
00036           .addVetoOnThisFinalState(zeefinder.clusteredPhotonsFinalState())
00037           .addVetoOnThisFinalState(zmmfinder.constituentsFinalState())
00038           .addVetoOnThisFinalState(zmmfinder.clusteredPhotonsFinalState());
00039       FastJets jetpro(jetinput, FastJets::KT, 0.7);
00040       addProjection(jetpro, "Jets");
00041 
00042       // properties of the pair momentum
00043       _h_ZZ_pT = bookHistogram1D("ZZ_pT", logBinEdges(100, 1.0, 0.5*sqrtS()));
00044       _h_ZZ_pT_peak = bookHistogram1D("ZZ_pT_peak", 25, 0.0, 25.0);
00045       _h_ZZ_eta = bookHistogram1D("ZZ_eta", 40, -7.0, 7.0);
00046       _h_ZZ_phi = bookHistogram1D("ZZ_phi", 25, 0.0, TWOPI);
00047       _h_ZZ_m = bookHistogram1D("ZZ_m", logBinEdges(100, 150.0, 180.0+0.25*sqrtS()));
00048 
00049       // correlations between the ZZ
00050       _h_ZZ_dphi = bookHistogram1D("ZZ_dphi", 25, 0.0, PI);  /// @todo non-linear?
00051       _h_ZZ_deta = bookHistogram1D("ZZ_deta", 25, -7.0, 7.0);
00052       _h_ZZ_dR = bookHistogram1D("ZZ_dR", 25, 0.5, 7.0);
00053       _h_ZZ_dpT = bookHistogram1D("ZZ_dpT", logBinEdges(100, 1.0, 0.5*sqrtS()));
00054       _h_ZZ_costheta_planes = bookHistogram1D("ZZ_costheta_planes", 25, -1.0, 1.0);
00055 
00056       /// @todo fuer WW: missing ET
00057 
00058       // properties of the Z bosons
00059       _h_Z_pT = bookHistogram1D("Z_pT", logBinEdges(100, 10.0, 0.25*sqrtS()));
00060       _h_Z_eta = bookHistogram1D("Z_eta", 70, -7.0, 7.0);
00061 
00062       // properties of the leptons
00063       _h_Zl_pT = bookHistogram1D("Zl_pT", logBinEdges(100, 30.0, 0.1
00064                                                       *sqrtS()));
00065       _h_Zl_eta = bookHistogram1D("Zl_eta", 40, -3.5, 3.5);
00066 
00067       // correlations between the opposite charge leptons
00068       _h_ZeZm_dphi = bookHistogram1D("ZeZm_dphi", 25, 0.0, PI);
00069       _h_ZeZm_deta = bookHistogram1D("ZeZm_deta", 25, -5.0, 5.0);
00070       _h_ZeZm_dR = bookHistogram1D("ZeZm_dR", 25, 0.5, 5.0);
00071       _h_ZeZm_m = bookHistogram1D("ZeZm_m", 100, 0.0, 300.0);
00072 
00073       // correlations with jets
00074       _h_ZZ_jet1_deta = bookHistogram1D("ZZ_jet1_deta", 70, -7.0, 7.0);
00075       _h_ZZ_jet1_dR = bookHistogram1D("ZZ_jet1_dR", 25, 1.5, 7.0);
00076       _h_Ze_jet1_dR = bookHistogram1D("Ze_jet1_dR", 25, 0.0, 7.0);
00077 
00078       // global stuff
00079       _h_HT = bookHistogram1D("HT", logBinEdges(100, 100.0, 0.5*sqrtS()));
00080 
00081       MC_JetAnalysis::init();
00082     }
00083 
00084 
00085 
00086     /// Do the analysis
00087     void analyze(const Event & e) {
00088       const double weight = e.weight();
00089 
00090       const ZFinder& zeefinder = applyProjection<ZFinder>(e, "ZeeFinder");
00091       if (zeefinder.particles().size()!=1) {
00092         vetoEvent;
00093       }
00094 
00095       const ZFinder& zmmfinder = applyProjection<ZFinder>(e, "ZmmFinder");
00096       if (zmmfinder.particles().size()!=1) {
00097         vetoEvent;
00098       }
00099 
00100       FourMomentum zee(zeefinder.particles()[0].momentum());
00101       FourMomentum zmm(zmmfinder.particles()[0].momentum());
00102       FourMomentum zz(zee+zmm);
00103       // find leptons
00104       FourMomentum ep(0.0,0.0,0.0,0.0), em(0.0,0.0,0.0,0.0);
00105       if (PID::threeCharge(zeefinder.constituentsFinalState().particles()[0])>0.0) {
00106         ep=zeefinder.constituentsFinalState().particles()[0].momentum();
00107         em=zeefinder.constituentsFinalState().particles()[1].momentum();
00108       }
00109       else {
00110         ep=zeefinder.constituentsFinalState().particles()[1].momentum();
00111         em=zeefinder.constituentsFinalState().particles()[0].momentum();
00112       }
00113       FourMomentum mp(0.0,0.0,0.0,0.0), mm(0.0,0.0,0.0,0.0);
00114       if (PID::threeCharge(zmmfinder.constituentsFinalState().particles()[0])>0.0) {
00115         mp=zmmfinder.constituentsFinalState().particles()[0].momentum();
00116         mm=zmmfinder.constituentsFinalState().particles()[1].momentum();
00117       }
00118       else {
00119         mp=zmmfinder.constituentsFinalState().particles()[1].momentum();
00120         mm=zmmfinder.constituentsFinalState().particles()[0].momentum();
00121       }
00122 
00123       _h_ZZ_pT->fill(zz.pT(),weight);
00124       _h_ZZ_pT_peak->fill(zz.pT(),weight);
00125       _h_ZZ_eta->fill(zz.eta(),weight);
00126       _h_ZZ_phi->fill(zz.azimuthalAngle(),weight);
00127       double mzz2=zz.mass2();
00128       if (mzz2>0.0) _h_ZZ_m->fill(sqrt(mzz2), weight);
00129 
00130       _h_ZZ_dphi->fill(mapAngle0ToPi(zee.phi()-zmm.phi()), weight);
00131       _h_ZZ_deta->fill(zee.eta()-zmm.eta(), weight);
00132       _h_ZZ_dR->fill(deltaR(zee,zmm), weight);
00133       _h_ZZ_dpT->fill(fabs(zee.pT()-zmm.pT()), weight);
00134 
00135       Vector3 crossZee = ep.vector3().cross(em.vector3());
00136       Vector3 crossZmm = mp.vector3().cross(mm.vector3());
00137       double costheta = crossZee.dot(crossZmm)/crossZee.mod()/crossZmm.mod();
00138       _h_ZZ_costheta_planes->fill(costheta, weight);
00139 
00140       _h_Z_pT->fill(zee.pT(),weight);
00141       _h_Z_pT->fill(zmm.pT(),weight);
00142       _h_Z_eta->fill(zee.eta(),weight);
00143       _h_Z_eta->fill(zmm.eta(),weight);
00144 
00145       _h_Zl_pT->fill(ep.pT(), weight);
00146       _h_Zl_pT->fill(em.pT(), weight);
00147       _h_Zl_pT->fill(mp.pT(), weight);
00148       _h_Zl_pT->fill(mm.pT(), weight);
00149       _h_Zl_eta->fill(ep.eta(), weight);
00150       _h_Zl_eta->fill(em.eta(), weight);
00151       _h_Zl_eta->fill(mp.eta(), weight);
00152       _h_Zl_eta->fill(mm.eta(), weight);
00153 
00154       _h_ZeZm_dphi->fill(mapAngle0ToPi(ep.phi()-mm.phi()), weight);
00155       _h_ZeZm_deta->fill(ep.eta()-mm.eta(), weight);
00156       _h_ZeZm_dR->fill(deltaR(ep,mm), weight);
00157       double m2=FourMomentum(ep+mm).mass2();
00158       if (m2 < 0) m2 = 0.0;
00159       _h_ZeZm_m->fill(sqrt(m2), weight);
00160 
00161       const FastJets& jetpro = applyProjection<FastJets>(e, "Jets");
00162       const Jets& jets = jetpro.jetsByPt(20.0*GeV);
00163       if (jets.size() > 0) {
00164         _h_ZZ_jet1_deta->fill(zz.eta()-jets[0].momentum().eta(), weight);
00165         _h_ZZ_jet1_dR->fill(deltaR(zz, jets[0].momentum()), weight);
00166         _h_Ze_jet1_dR->fill(deltaR(ep, jets[0].momentum()), weight);
00167       }
00168 
00169       double HT=ep.pT()+em.pT()+mp.pT()+mm.pT();
00170       foreach (const Jet& jet, jets) {
00171         HT+=jet.momentum().pT();
00172       }
00173       foreach (const Particle& p, zeefinder.clusteredPhotonsFinalState().particles()) {
00174         HT+=p.momentum().pT();
00175       }
00176       foreach (const Particle& p, zmmfinder.clusteredPhotonsFinalState().particles()) {
00177         HT+=p.momentum().pT();
00178       }
00179       if (HT>0.0) _h_HT->fill(HT, weight);
00180 
00181       MC_JetAnalysis::analyze(e);
00182     }
00183 
00184 
00185     /// Finalize
00186     void finalize() {
00187       double norm=crossSection()/sumOfWeights();
00188       scale(_h_ZZ_pT, norm);
00189       scale(_h_ZZ_pT_peak, norm);
00190       scale(_h_ZZ_eta, norm);
00191       scale(_h_ZZ_phi, norm);
00192       scale(_h_ZZ_m, norm);
00193       scale(_h_ZZ_dphi, norm);
00194       scale(_h_ZZ_deta, norm);
00195       scale(_h_ZZ_dR, norm);
00196       scale(_h_ZZ_dpT, norm);
00197       scale(_h_ZZ_costheta_planes, norm);
00198       scale(_h_Z_pT, norm);
00199       scale(_h_Z_eta, norm);
00200       scale(_h_Zl_pT, norm);
00201       scale(_h_Zl_eta, norm);
00202       scale(_h_ZeZm_dphi, norm);
00203       scale(_h_ZeZm_deta, norm);
00204       scale(_h_ZeZm_dR, norm);
00205       scale(_h_ZeZm_m, norm);
00206       scale(_h_ZZ_jet1_deta, norm);
00207       scale(_h_ZZ_jet1_dR, norm);
00208       scale(_h_Ze_jet1_dR, norm);
00209       scale(_h_HT, norm);
00210 
00211       MC_JetAnalysis::finalize();
00212     }
00213 
00214     //@}
00215 
00216 
00217   private:
00218 
00219     /// @name Histograms
00220     //@{
00221     AIDA::IHistogram1D * _h_ZZ_pT;
00222     AIDA::IHistogram1D * _h_ZZ_pT_peak;
00223     AIDA::IHistogram1D * _h_ZZ_eta;
00224     AIDA::IHistogram1D * _h_ZZ_phi;
00225     AIDA::IHistogram1D * _h_ZZ_m;
00226     AIDA::IHistogram1D * _h_ZZ_dphi;
00227     AIDA::IHistogram1D * _h_ZZ_deta;
00228     AIDA::IHistogram1D * _h_ZZ_dR;
00229     AIDA::IHistogram1D * _h_ZZ_dpT;
00230     AIDA::IHistogram1D * _h_ZZ_costheta_planes;
00231     AIDA::IHistogram1D * _h_Z_pT;
00232     AIDA::IHistogram1D * _h_Z_eta;
00233     AIDA::IHistogram1D * _h_Zl_pT;
00234     AIDA::IHistogram1D * _h_Zl_eta;
00235     AIDA::IHistogram1D * _h_ZeZm_dphi;
00236     AIDA::IHistogram1D * _h_ZeZm_deta;
00237     AIDA::IHistogram1D * _h_ZeZm_dR;
00238     AIDA::IHistogram1D * _h_ZeZm_m;
00239     AIDA::IHistogram1D * _h_ZZ_jet1_deta;
00240     AIDA::IHistogram1D * _h_ZZ_jet1_dR;
00241     AIDA::IHistogram1D * _h_Ze_jet1_dR;
00242     AIDA::IHistogram1D * _h_HT;
00243     //@}
00244 
00245   };
00246 
00247 
00248 
00249   // This global object acts as a hook for the plugin system
00250   AnalysisBuilder<MC_ZZJETS> plugin_MC_ZZJETS;
00251 
00252 }