[Rivet-svn] r1880 - in trunk: bin src/Analyses

[blackhole at projects.hepforge.org]

Author: buckley
Date: Tue Oct  6 18:07:21 2009
New Revision: 1880

Log:
Bug fixes and isolated photons in SUSY analysis

Modified:
   trunk/bin/rivet
   trunk/src/Analyses/MC_LHC_SUSY.cc
   trunk/src/Analyses/OPAL_2004_S6132243.cc

Modified: trunk/bin/rivet
==============================================================================
--- trunk/bin/rivet	Tue Oct  6 16:39:28 2009	(r1879)
+++ trunk/bin/rivet	Tue Oct  6 18:07:21 2009	(r1880)
@@ -142,8 +142,8 @@
             level = rivet.Log.DEBUG
         elif LEVEL == "INFO":
             level = rivet.Log.INFO
-        elif LEVEL == "WARNING":
-            level = rivet.Log.WARNING
+        elif LEVEL == "WARNING" or LEVEL == "WARN":
+            level = rivet.Log.WARN
         elif LEVEL == "ERROR":
             level = rivet.Log.ERROR
         else:

Modified: trunk/src/Analyses/MC_LHC_SUSY.cc
==============================================================================
--- trunk/src/Analyses/MC_LHC_SUSY.cc	Tue Oct  6 16:39:28 2009	(r1879)
+++ trunk/src/Analyses/MC_LHC_SUSY.cc	Tue Oct  6 18:07:21 2009	(r1880)
@@ -82,7 +82,10 @@
       _hist_eta_gamma = bookHistogram1D("eta-gamma", 50, -4, 4);
       _hist_pt_gamma  = bookHistogram1D("pt-gamma", 100, 0.0, 500);
 
-      /// @todo Isolated photons
+      _hist_n_gammaiso   = bookHistogram1D("n-gamma-iso", 11, -0.5, 10.5);
+      _hist_phi_gammaiso = bookHistogram1D("phi-gamma-iso", 50, -PI, PI);
+      _hist_eta_gammaiso = bookHistogram1D("eta-gamma-iso", 50, -4, 4);
+      _hist_pt_gammaiso  = bookHistogram1D("pt-gamma-iso", 100, 0.0, 500);
 
       _hist_met = bookHistogram1D("Etmiss", 100, 0.0, 1500);
 
@@ -109,7 +112,7 @@
       _hist_n_trk->fill(tracks.size(), weight);
       foreach (const Particle& t, tracks.particles()) {
         const FourMomentum& p = t.momentum();
-        _hist_phi_trk->fill(p.phi(), weight);
+        _hist_phi_trk->fill(mapAngleMPiToPi(p.phi()), weight);
         _hist_eta_trk->fill(p.eta(), weight);
         _hist_pt_trk->fill(p.pT()/GeV, weight);
       }
@@ -121,7 +124,7 @@
       _hist_n_jet->fill(jets.size(), weight);
       foreach (const Jet& j, jets) {
         const FourMomentum& pj = j.momentum();
-        _hist_phi_jet->fill(pj.phi(), weight);
+        _hist_phi_jet->fill(mapAngleMPiToPi(pj.phi()), weight);
         _hist_eta_jet->fill(pj.eta(), weight);
         _hist_pt_jet->fill(pj.pT()/GeV, weight);
       }
@@ -133,7 +136,7 @@
       _hist_n_e->fill(efs.size(), weight);
       foreach (const Particle& e, efs.particles()) {
         const FourMomentum& p = e.momentum();
-        _hist_phi_e->fill(p.phi(), weight);
+        _hist_phi_e->fill(mapAngleMPiToPi(p.phi()), weight);
         _hist_eta_e->fill(p.eta(), weight);
         _hist_pt_e->fill(p.pT()/GeV, weight);
       }
@@ -142,25 +145,43 @@
 
       // Fill final state muon/antimuon histos
       const FinalState& mufs = applyProjection<FinalState>(evt, "Muons");
-      _hist_n_e->fill(mufs.size(), weight);
+      _hist_n_mu->fill(mufs.size(), weight);
       foreach (const Particle& mu, efs.particles()) {
         const FourMomentum& p = mu.momentum();
-        _hist_phi_mu->fill(p.phi(), weight);
+        _hist_phi_mu->fill(mapAngleMPiToPi(p.phi()), weight);
         _hist_eta_mu->fill(p.eta(), weight);
         _hist_pt_mu->fill(p.pT()/GeV, weight);
       }
 
       // Fill final state non-isolated photon histos
       const FinalState& allphotonfs = applyProjection<FinalState>(evt, "AllPhotons");
-      _hist_n_e->fill(allphotonfs.size(), weight);
+      _hist_n_gamma->fill(allphotonfs.size(), weight);
+      ParticleVector isolatedphotons;
       foreach (const Particle& ph, allphotonfs.particles()) {
         const FourMomentum& p = ph.momentum();
-        _hist_phi_mu->fill(p.phi(), weight);
-        _hist_eta_mu->fill(p.eta(), weight);
-        _hist_pt_mu->fill(p.pT()/GeV, weight);
+        _hist_phi_gamma->fill(mapAngleMPiToPi(p.phi()), weight);
+        _hist_eta_gamma->fill(p.eta(), weight);
+        _hist_pt_gamma->fill(p.pT()/GeV, weight);
+        // Select isolated photons
+        bool isolated = true;
+        foreach (const Jet& j, jets) {
+          if (deltaR(j.momentum(), p) < 0.2) {
+            isolated = false;
+            break;
+          }
+        }
+        if (isolated) isolatedphotons += ph;
       }
 
-      /// @todo Isolated photons
+
+      // Fill final state isolated photon histos
+      _hist_n_gammaiso->fill(isolatedphotons.size(), weight);
+      foreach (const Particle& ph_iso, isolatedphotons) {
+        const FourMomentum& p = ph_iso.momentum();
+        _hist_phi_gammaiso->fill(mapAngleMPiToPi(p.phi()), weight);
+        _hist_eta_gammaiso->fill(p.eta(), weight);
+        _hist_pt_gammaiso->fill(p.pT()/GeV, weight);
+      }
 
       // Calculate and fill missing Et histos
       const TotalVisibleMomentum& met = applyProjection<TotalVisibleMomentum>(evt, "MET");
@@ -225,7 +246,7 @@
     AIDA::IHistogram1D *_hist_n_e, *_hist_phi_e, *_hist_eta_e, *_hist_pt_e;
     AIDA::IHistogram1D *_hist_n_mu, *_hist_phi_mu, *_hist_eta_mu, *_hist_pt_mu;
     AIDA::IHistogram1D *_hist_n_gamma, *_hist_phi_gamma, *_hist_eta_gamma, *_hist_pt_gamma;
-    /// @todo Isolated photons
+    AIDA::IHistogram1D *_hist_n_gammaiso, *_hist_phi_gammaiso, *_hist_eta_gammaiso, *_hist_pt_gammaiso;
     AIDA::IHistogram1D *_hist_met;
     AIDA::IHistogram1D *_hist_mll_ossf_ee, *_hist_mll_ossf_mumu, *_hist_mll_osof_emu;
     

Modified: trunk/src/Analyses/OPAL_2004_S6132243.cc
==============================================================================
--- trunk/src/Analyses/OPAL_2004_S6132243.cc	Tue Oct  6 16:39:28 2009	(r1879)
+++ trunk/src/Analyses/OPAL_2004_S6132243.cc	Tue Oct  6 18:07:21 2009	(r1880)
@@ -75,7 +75,7 @@
 
 
     void analyze(const Event& event) { 
-      const FinalState& cfs = applyProjection<FinalState>(event, "CFS");
+      const FinalState& cfs = applyProjection<FinalState>(event, "FS");
       // Even if we only generate hadronic events, we still need a cut on numCharged >= 2.
       if (cfs.size() < 2) vetoEvent;
 
@@ -184,7 +184,7 @@
         normalize(_histDParam[isqrts]);
         normalize(_histY23Durham[isqrts]);
         //
-        scale(_hist1MinusTMom[isqrts], 1.0/_sumPassedWeights); 
+        scale(_hist1MinusTMom[isqrts], 1.0/_sumPassedWeights);
         scale(_histTMajorMom[isqrts], 1.0/_sumPassedWeights); 
         scale(_histTMinorMom[isqrts], 1.0/_sumPassedWeights); 
         scale(_histOblatenessMom[isqrts], 1.0/_sumPassedWeights);