ATLAS_2010_S8918562.cc

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 
00003 #include "Rivet/Analysis.hh"
00004 #include "Rivet/RivetAIDA.hh"
00005 #include "Rivet/Projections/ChargedFinalState.hh"
00006 #include "Rivet/Tools/Logging.hh"
00007 
00008 namespace Rivet {
00009 
00010 
00011   /// Rivet analysis class for ATLAS_2010_S8918562 dataset
00012   class ATLAS_2010_S8918562 : public Analysis {
00013   public:
00014 
00015     /// Helper for collectively filling Nch, pT, eta, and pT vs. Nch histograms
00016     void fillPtEtaNch(const ChargedFinalState& cfs, int nchcut, double weight,
00017                       AIDA::IHistogram1D* h_nch, AIDA::IHistogram1D* h_pt,
00018                       AIDA::IHistogram1D* h_eta, AIDA::IProfile1D* h_ptnch = 0) {
00019 
00020       // Get number of particles and skip if event fails cut
00021       const int nch = cfs.size();
00022       if (nch < nchcut) return;
00023 
00024       // Fill nch
00025       h_nch->fill(nch, weight);
00026       // Loop over particles, fill pT, eta and ptnch
00027       foreach (const Particle& p, cfs.particles()) {
00028         const double pt = p.momentum().pT();
00029         h_pt->fill(pt/GeV, weight/pt);
00030         h_eta->fill(p.momentum().eta(), weight);
00031         if (h_ptnch != 0) h_ptnch->fill(nch, pt/GeV, weight);
00032       }
00033     }
00034 
00035 
00036     /// Default constructor
00037     ATLAS_2010_S8918562() : Analysis("ATLAS_2010_S8918562") {
00038       _sumW_pt100_nch2 = 0;
00039       _sumW_pt100_nch20 = 0;
00040       _sumW_pt500_nch1 = 0;
00041       _sumW_pt500_nch6 = 0;
00042       _sumW_pt2500_nch1 = 0;
00043     }
00044 
00045 
00046     /// Initialization, called once before running
00047     void init() {
00048       // Projections
00049       const ChargedFinalState cfs100(-2.5, 2.5, 100.0*MeV);
00050       addProjection(cfs100, "CFS100");
00051       const ChargedFinalState cfs500(-2.5, 2.5, 500.0*MeV);
00052       addProjection(cfs500, "CFS500");
00053       const ChargedFinalState cfs2500(-2.5, 2.5, 2500.0*MeV);
00054       addProjection(cfs2500, "CFS2500");
00055 
00056       // Book histograms
00057       if (fuzzyEquals(sqrtS()/GeV, 900)) {
00058         _hist_pt100_nch2_nch = bookHistogram1D(18, 1, 1);
00059         _hist_pt100_nch2_pt = bookHistogram1D(11, 1, 1);
00060         _hist_pt100_nch2_eta = bookHistogram1D(4, 1, 1);
00061         _hist_pt100_nch2_ptnch = bookProfile1D(24, 1, 1);
00062 
00063         _hist_pt100_nch20_nch = bookHistogram1D(34, 1, 1);
00064         _hist_pt100_nch20_pt = bookHistogram1D(30, 1, 1);
00065         _hist_pt100_nch20_eta = bookHistogram1D(26, 1, 1);
00066 
00067         _hist_pt500_nch1_nch = bookHistogram1D(15, 1, 1);
00068         _hist_pt500_nch1_pt = bookHistogram1D(8, 1, 1);
00069         _hist_pt500_nch1_eta = bookHistogram1D(1, 1, 1);
00070         _hist_pt500_nch1_ptnch = bookProfile1D(22, 1, 1);
00071 
00072         _hist_pt500_nch6_nch = bookHistogram1D(20, 1, 1);
00073         _hist_pt500_nch6_pt = bookHistogram1D(13, 1, 1);
00074         _hist_pt500_nch6_eta = bookHistogram1D(6, 1, 1);
00075 
00076         _hist_pt2500_nch1_nch = bookHistogram1D(36, 1, 1);
00077         _hist_pt2500_nch1_pt = bookHistogram1D(32, 1, 1);
00078         _hist_pt2500_nch1_eta = bookHistogram1D(28, 1, 1);
00079         _hist_pt2500_nch1_ptnch = bookProfile1D(38, 1, 1);
00080 
00081       } else if (fuzzyEquals(sqrtS()/GeV, 2360)) {
00082 
00083         _hist_pt500_nch1_nch = bookHistogram1D(16, 1, 1);
00084         _hist_pt500_nch1_pt = bookHistogram1D(9, 1, 1);
00085         _hist_pt500_nch1_eta = bookHistogram1D(2, 1, 1);
00086         // This one histogram might be called while unbooked, so ensure its pointer is null!
00087         _hist_pt500_nch1_ptnch = 0;
00088 
00089       } else if (fuzzyEquals(sqrtS()/GeV, 7000)) {
00090 
00091         _hist_pt100_nch2_nch = bookHistogram1D(19, 1, 1);
00092         _hist_pt100_nch2_pt = bookHistogram1D(12, 1, 1);
00093         _hist_pt100_nch2_eta = bookHistogram1D(5, 1, 1);
00094         _hist_pt100_nch2_ptnch = bookProfile1D(25, 1, 1);
00095 
00096         _hist_pt100_nch20_nch = bookHistogram1D(35, 1, 1);
00097         _hist_pt100_nch20_pt = bookHistogram1D(31, 1, 1);
00098         _hist_pt100_nch20_eta = bookHistogram1D(27, 1, 1);
00099 
00100         _hist_pt500_nch1_nch = bookHistogram1D(17, 1, 1);
00101         _hist_pt500_nch1_pt = bookHistogram1D(10, 1, 1);
00102         _hist_pt500_nch1_eta = bookHistogram1D(3, 1, 1);
00103         _hist_pt500_nch1_ptnch = bookProfile1D(23, 1, 1);
00104 
00105         _hist_pt500_nch6_nch = bookHistogram1D(21, 1, 1);
00106         _hist_pt500_nch6_pt = bookHistogram1D(14, 1, 1);
00107         _hist_pt500_nch6_eta = bookHistogram1D(7, 1, 1);
00108 
00109         _hist_pt2500_nch1_nch = bookHistogram1D(37, 1, 1);
00110         _hist_pt2500_nch1_pt = bookHistogram1D(33, 1, 1);
00111         _hist_pt2500_nch1_eta = bookHistogram1D(29, 1, 1);
00112         _hist_pt2500_nch1_ptnch = bookProfile1D(39, 1, 1);
00113 
00114       } else {
00115         throw LogicError("The ATLAS_2010_S8918562 analysis is only valid for sqrt(s) = 900, 2360 and 7000 GeV!");
00116       }
00117 
00118     }
00119 
00120 
00121 
00122     void analyze(const Event& event) {
00123       const double weight = event.weight();
00124 
00125       // 100 GeV final states
00126       if (!fuzzyEquals(sqrtS()/GeV, 2360)) {
00127         const ChargedFinalState& cfs100 = applyProjection<ChargedFinalState>(event, "CFS100");
00128         // nch>=2
00129         if (cfs100.size() >= 2) _sumW_pt100_nch2 += weight;
00130         fillPtEtaNch(cfs100, 2, weight, _hist_pt100_nch2_nch, _hist_pt100_nch2_pt, _hist_pt100_nch2_eta, _hist_pt100_nch2_ptnch);
00131         // nch>=20
00132         if (cfs100.size() >= 20) _sumW_pt100_nch20 += weight;
00133         fillPtEtaNch(cfs100, 20, weight, _hist_pt100_nch20_nch, _hist_pt100_nch20_pt, _hist_pt100_nch20_eta);
00134       }
00135 
00136       // 500 GeV final states
00137       const ChargedFinalState& cfs500 = applyProjection<ChargedFinalState>(event, "CFS500");
00138       // nch>=1
00139       if (cfs500.size() >= 1) _sumW_pt500_nch1 += weight;
00140       fillPtEtaNch(cfs500, 1, weight, _hist_pt500_nch1_nch, _hist_pt500_nch1_pt, _hist_pt500_nch1_eta, _hist_pt500_nch1_ptnch);
00141       // nch>=6
00142       if (!fuzzyEquals(sqrtS()/GeV, 2360)) {
00143         if (cfs500.size() >= 6) _sumW_pt500_nch6 += weight;
00144         fillPtEtaNch(cfs500, 6, weight, _hist_pt500_nch6_nch, _hist_pt500_nch6_pt, _hist_pt500_nch6_eta);
00145       }
00146 
00147       // 2500 GeV final states
00148       if (!fuzzyEquals(sqrtS()/GeV, 2360)) {
00149         const ChargedFinalState& cfs2500 = applyProjection<ChargedFinalState>(event, "CFS2500");
00150         // nch>=1
00151         if (cfs2500.size() >= 1) _sumW_pt2500_nch1 += weight;
00152         fillPtEtaNch(cfs2500, 1, weight, _hist_pt2500_nch1_nch, _hist_pt2500_nch1_pt, _hist_pt2500_nch1_eta, _hist_pt2500_nch1_ptnch);
00153       }
00154 
00155     }
00156 
00157 
00158 
00159     void finalize() {
00160 
00161       if (!fuzzyEquals(sqrtS()/GeV, 2360)) {
00162         if (_sumW_pt100_nch2 > 0) {
00163           _hist_pt100_nch2_nch->scale(1.0/_sumW_pt100_nch2);
00164           _hist_pt100_nch2_pt->scale(1.0/_sumW_pt100_nch2/TWOPI/5);
00165           _hist_pt100_nch2_eta->scale(1.0/_sumW_pt100_nch2);
00166         }
00167 
00168         if (_sumW_pt100_nch20 > 0) {
00169           _hist_pt100_nch20_nch->scale(1.0/_sumW_pt100_nch20);
00170           _hist_pt100_nch20_pt->scale(1.0/_sumW_pt100_nch20/TWOPI/5);
00171           _hist_pt100_nch20_eta->scale(1.0/_sumW_pt100_nch20);
00172         }
00173 
00174         if (_sumW_pt500_nch6 > 0) {
00175           _hist_pt500_nch6_nch->scale(1.0/_sumW_pt500_nch6);
00176           _hist_pt500_nch6_pt->scale(1.0/_sumW_pt500_nch6/TWOPI/5);
00177           _hist_pt500_nch6_eta->scale(1.0/_sumW_pt500_nch6);
00178         }
00179 
00180         if (_sumW_pt2500_nch1 > 0) {
00181           _hist_pt2500_nch1_nch->scale(1.0/_sumW_pt2500_nch1);
00182           _hist_pt2500_nch1_pt->scale(1.0/_sumW_pt2500_nch1/TWOPI/5);
00183           _hist_pt2500_nch1_eta->scale(1.0/_sumW_pt2500_nch1);
00184         }
00185       }
00186 
00187       if (_sumW_pt500_nch1 > 0) {
00188         _hist_pt500_nch1_nch->scale(1.0/_sumW_pt500_nch1);
00189         _hist_pt500_nch1_pt->scale(1.0/_sumW_pt500_nch1/TWOPI/5);
00190         _hist_pt500_nch1_eta->scale(1.0/_sumW_pt500_nch1);
00191       }
00192     }
00193 
00194 
00195   private:
00196 
00197     double _sumW_pt100_nch2, _sumW_pt100_nch20, _sumW_pt500_nch1,
00198       _sumW_pt500_nch6, _sumW_pt2500_nch1;
00199 
00200     AIDA::IHistogram1D *_hist_pt100_nch2_nch,
00201       *_hist_pt100_nch20_nch, *_hist_pt500_nch1_nch,
00202       *_hist_pt500_nch6_nch, *_hist_pt2500_nch1_nch;
00203 
00204     AIDA::IHistogram1D *_hist_pt100_nch2_pt,
00205       *_hist_pt100_nch20_pt, *_hist_pt500_nch1_pt,
00206       *_hist_pt500_nch6_pt, *_hist_pt2500_nch1_pt;
00207 
00208     AIDA::IHistogram1D *_hist_pt100_nch2_eta,
00209       *_hist_pt100_nch20_eta, *_hist_pt500_nch1_eta,
00210       *_hist_pt500_nch6_eta, *_hist_pt2500_nch1_eta;
00211 
00212     AIDA::IProfile1D *_hist_pt100_nch2_ptnch,
00213       *_hist_pt500_nch1_ptnch, *_hist_pt2500_nch1_ptnch;
00214 
00215   };
00216 
00217 
00218 
00219   // This global object acts as a hook for the plugin system
00220   AnalysisBuilder<ATLAS_2010_S8918562> plugin_ATLAS_2010_S8918562;
00221 
00222 }