1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/FinalState.hh"
  4#include "Rivet/Projections/ChargedFinalState.hh"
  5#include "Rivet/Projections/VisibleFinalState.hh"
  6#include "Rivet/Projections/VetoedFinalState.hh"
  7#include "Rivet/Projections/IdentifiedFinalState.hh"
  8#include "Rivet/Projections/FastJets.hh"
  9#include "Rivet/Tools/RivetMT2.hh"
 10
 11namespace Rivet {
 12
 13
 14  class ATLAS_2012_CONF_2012_103 : public Analysis {
 15  public:
 16
 17    /// Constructor
 18    ATLAS_2012_CONF_2012_103()
 19      : Analysis("ATLAS_2012_CONF_2012_103")
 20    {    }
 21
 22
 23    /// @name Analysis methods
 24    /// @{
 25
 26    /// Book histograms and initialise projections before the run
 27    void init() {
 28
 29      // projection to find the electrons
 30      IdentifiedFinalState elecs(Cuts::abseta < 2.47 && Cuts::pT > 20*GeV);
 31      elecs.acceptIdPair(PID::ELECTRON);
 32      declare(elecs, "elecs");
 33
 34      // projection to find the muons
 35      IdentifiedFinalState muons(Cuts::abseta < 2.4 && Cuts::pT > 10*GeV);
 36      muons.acceptIdPair(PID::MUON);
 37      declare(muons, "muons");
 38
 39      // for pTmiss
 40      declare(VisibleFinalState(Cuts::abseta < 4.9), "vfs");
 41
 42      VetoedFinalState vfs;
 43      vfs.addVetoPairId(PID::MUON);
 44
 45      /// Jet finder
 46      declare(FastJets(vfs, JetAlg::ANTIKT, 0.4), "AntiKtJets04");
 47
 48      /// Book histograms
 49      book(_etmiss_HT_7j55 ,"etmiss_HT_7j55", 8, 0., 16.);
 50      book(_etmiss_HT_8j55 ,"etmiss_HT_8j55", 8, 0., 16.);
 51      book(_etmiss_HT_9j55 ,"etmiss_HT_9j55", 8, 0., 16.);
 52      book(_etmiss_HT_6j80 ,"etmiss_HT_6j80", 8, 0., 16.);
 53      book(_etmiss_HT_7j80 ,"etmiss_HT_7j80", 8, 0., 16.);
 54      book(_etmiss_HT_8j80 ,"etmiss_HT_8j80", 8, 0., 16.);
 55
 56      book(_hist_njet55 ,"hist_njet55", 4, 5.5, 9.5);
 57      book(_hist_njet80 ,"hist_njet80", 4, 4.5, 8.5);
 58
 59      book(_count_7j55 ,"count_7j55", 1, 0., 1.);
 60      book(_count_8j55 ,"count_8j55", 1, 0., 1.);
 61      book(_count_9j55 ,"count_9j55", 1, 0., 1.);
 62      book(_count_6j80 ,"count_6j80", 1, 0., 1.);
 63      book(_count_7j80 ,"count_7j80", 1, 0., 1.);
 64      book(_count_8j80 ,"count_8j80", 1, 0., 1.);
 65
 66    }
 67
 68
 69    /// Perform the per-event analysis
 70    void analyze(const Event& event) {
 71
 72      // get the jet candidates
 73      Jets cand_jets = apply<FastJets>(event, "AntiKtJets04").jetsByPt(Cuts::pT > 20*GeV && Cuts::abseta < 2.8);
 74
 75      // candidate muons
 76      Particles cand_mu = apply<IdentifiedFinalState>(event, "muons").particlesByPt();
 77
 78      // candidate electrons
 79      Particles cand_e  = apply<IdentifiedFinalState>(event, "elecs").particlesByPt();
 80
 81      // resolve jet/lepton ambiguity
 82      Jets recon_jets;
 83      for ( const Jet& jet : cand_jets ) {
 84        // candidates after |eta| < 2.8
 85        if ( fabs( jet.eta() ) >= 2.8 ) continue;
 86        bool away_from_e = true;
 87        for ( const Particle & e : cand_e ) {
 88          if ( deltaR(e.momentum(),jet.momentum()) <= 0.2 ) {
 89            away_from_e = false;
 90            break;
 91          }
 92        }
 93        if ( away_from_e ) recon_jets.push_back( jet );
 94      }
 95
 96      // only keep electrons more than R=0.4 from jets
 97      Particles recon_e;
 98      for ( const Particle & e : cand_e ) {
 99        bool away = true;
100        for ( const Jet& jet : recon_jets ) {
101          if ( deltaR(e.momentum(),jet.momentum()) < 0.4 ) {
102            away = false;
103            break;
104          }
105        }
106        if ( away )
107          recon_e.push_back( e );
108      }
109
110      // only keep muons more than R=0.4 from jets
111      Particles recon_mu;
112      for ( const Particle & mu : cand_mu ) {
113        bool away = true;
114        for ( const Jet& jet : recon_jets ) {
115          if ( deltaR(mu.momentum(),jet.momentum()) < 0.4 ) {
116            away = false;
117            break;
118          }
119        }
120        if ( away )
121          recon_mu.push_back( mu );
122      }
123
124      // pTmiss
125      Particles vfs_particles =
126        apply<VisibleFinalState>(event, "vfs").particles();
127      FourMomentum pTmiss;
128      for ( const Particle & p : vfs_particles ) {
129        pTmiss -= p.momentum();
130      }
131      double eTmiss = pTmiss.pT();
132
133      // now only use recon_jets, recon_mu, recon_e
134
135      // reject events with electrons and muons
136      if ( ! ( recon_mu.empty() && recon_e.empty() ) ) {
137        MSG_DEBUG("Charged leptons left after selection");
138        vetoEvent;
139      }
140
141      // calculate H_T
142      double HT=0;
143      for ( const Jet& jet : recon_jets ) {
144        if ( jet.pT() > 40 * GeV )
145          HT += jet.pT() ;
146      }
147
148      // number of jets
149      unsigned int njet55=0, njet80=0;
150      for (unsigned int ix=0;ix<recon_jets.size();++ix) {
151        if(recon_jets[ix].pT()>80.*GeV) ++njet80;
152        if(recon_jets[ix].pT()>55.*GeV) ++njet55;
153      }
154
155      double ratio = eTmiss/sqrt(HT);
156
157      if(ratio>4.) {
158        if(njet55>9) njet55 = 9;
159        if(njet80>8) njet80 = 8;
160        _hist_njet55->fill(njet55);
161        _hist_njet80->fill(njet80);
162        // 7j55
163        if(njet55>=7)
164          _count_7j55->fill( 0.5);
165        // 8j55
166        if(njet55>=8)
167          _count_8j55->fill( 0.5) ;
168        // 8j55
169        if(njet55==9)
170          _count_9j55->fill( 0.5) ;
171        // 6j80
172        if(njet80>=6)
173          _count_6j80->fill( 0.5) ;
174        // 7j80
175        if(njet80>=7)
176          _count_7j80->fill( 0.5) ;
177        // 8j80
178        if(njet80==8)
179          _count_8j80->fill( 0.5) ;
180      }
181
182      if(njet55>=7)
183        _etmiss_HT_7j55->fill( ratio);
184      // 8j55
185      if(njet55>=8)
186        _etmiss_HT_8j55->fill( ratio) ;
187      // 8j55
188      if(njet55>=9)
189        _etmiss_HT_9j55->fill( ratio) ;
190      // 6j80
191      if(njet80>=6)
192        _etmiss_HT_6j80->fill( ratio) ;
193      // 7j80
194      if(njet80>=7)
195        _etmiss_HT_7j80->fill( ratio) ;
196      // 8j80
197      if(njet80>=8)
198        _etmiss_HT_8j80->fill( ratio) ;
199
200    }
201
202    /// @}
203
204    void finalize() {
205      double norm = crossSection()/femtobarn*5.8/sumOfWeights();
206
207      scale(_etmiss_HT_7j55,2.*norm);
208      scale(_etmiss_HT_8j55,2.*norm);
209      scale(_etmiss_HT_9j55,2.*norm);
210      scale(_etmiss_HT_6j80,2.*norm);
211      scale(_etmiss_HT_7j80,2.*norm);
212      scale(_etmiss_HT_8j80,2.*norm);
213
214      scale(_hist_njet55,norm);
215      scale(_hist_njet80,norm);
216
217      scale(_count_7j55,norm);
218      scale(_count_8j55,norm);
219      scale(_count_9j55,norm);
220      scale(_count_6j80,norm);
221      scale(_count_7j80,norm);
222      scale(_count_8j80,norm);
223    }
224
225  private:
226
227    /// @name Histograms
228    /// @{
229    Histo1DPtr _etmiss_HT_7j55;
230    Histo1DPtr _etmiss_HT_8j55;
231    Histo1DPtr _etmiss_HT_9j55;
232    Histo1DPtr _etmiss_HT_6j80;
233    Histo1DPtr _etmiss_HT_7j80;
234    Histo1DPtr _etmiss_HT_8j80;
235
236    Histo1DPtr _hist_njet55;
237    Histo1DPtr _hist_njet80;
238
239    Histo1DPtr _count_7j55;
240    Histo1DPtr _count_8j55;
241    Histo1DPtr _count_9j55;
242    Histo1DPtr _count_6j80;
243    Histo1DPtr _count_7j80;
244    Histo1DPtr _count_8j80;
245    /// @}
246
247  };
248
249  RIVET_DECLARE_PLUGIN(ATLAS_2012_CONF_2012_103);
250
251}