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/IdentifiedFinalState.hh"
  7#include "Rivet/Projections/FastJets.hh"
  8#include "Rivet/Projections/VetoedFinalState.hh"
  9
 10namespace Rivet {
 11
 12
 13  /// Two-lepton supersymmetry search
 14  class ATLAS_2011_I894578 : public Analysis {
 15  public:
 16
 17    /// Constructor
 18    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2011_I894578);
 19
 20
 21    /// @name Analysis methods
 22    /// @{
 23
 24    /// Book histograms and initialise projections before the run
 25    void init() {
 26
 27      // projection to find the electrons
 28      IdentifiedFinalState elecs(Cuts::abseta < 2.47 && Cuts::pT > 20*GeV);
 29      elecs.acceptIdPair(PID::ELECTRON);
 30      declare(elecs, "elecs");
 31
 32
 33      // veto region electrons
 34      Cut vetocut = Cuts::absetaIn(1.37, 1.52);
 35      IdentifiedFinalState veto_elecs(vetocut && Cuts::pT > 10*GeV);
 36      veto_elecs.acceptIdPair(PID::ELECTRON);
 37      declare(veto_elecs, "veto_elecs");
 38
 39
 40      // projection to find the muons
 41      IdentifiedFinalState muons(Cuts::abseta < 2.4 && Cuts::pT > 20*GeV);
 42      muons.acceptIdPair(PID::MUON);
 43      declare(muons, "muons");
 44
 45
 46      // jet finder
 47      VetoedFinalState vfs;
 48      vfs.addVetoPairId(PID::MUON);
 49      declare(FastJets(vfs, JetAlg::ANTIKT, 0.4), "AntiKtJets04");
 50
 51
 52      // all tracks (to do deltaR with leptons)
 53      declare(ChargedFinalState(Cuts::abseta < 3 && Cuts::pT > 0.5*GeV), "cfs");
 54
 55
 56      // for pTmiss
 57      declare(VisibleFinalState(Cuts::abseta < 4.9),"vfs");
 58
 59
 60      /// book histograms
 61      book(_count_OS_e_mu ,"count_OS_e+-mu-+", 1, 0., 1.);
 62      book(_count_OS_e_e ,"count_OS_e+e-", 1, 0., 1.);
 63      book(_count_OS_mu_mu ,"count_OS_mu+mu-", 1, 0., 1.);
 64      book(_count_SS_e_mu ,"count_SS_e+-mu+-", 1, 0., 1.);
 65      book(_count_SS_e_e ,"count_SS_e+-e+-", 1, 0., 1.);
 66      book(_count_SS_mu_mu ,"count_SS_mu+-mu+-", 1, 0., 1.);
 67
 68      book(_hist_eTmiss_OS  ,"Et_miss_OS", 20, 0., 400.);
 69      book(_hist_eTmiss_SS  ,"Et_miss_SS", 20, 0., 400.);
 70
 71    }
 72
 73
 74
 75
 76    /// Perform the per-event analysis
 77    void analyze(const Event& event) {
 78
 79
 80      Particles veto_e
 81        = apply<IdentifiedFinalState>(event, "veto_elecs").particles();
 82      if ( ! veto_e.empty() ) {
 83        MSG_DEBUG("electrons in veto region");
 84        vetoEvent;
 85      }
 86
 87      Jets cand_jets = apply<FastJets>(event, "AntiKtJets04").jetsByPt(Cuts::pT>20*GeV && Cuts::abseta < 2.5);
 88
 89      Particles cand_e =
 90        apply<IdentifiedFinalState>(event, "elecs").particlesByPt();
 91
 92      // charged particle for isolation
 93      Particles chg_tracks =
 94        apply<ChargedFinalState>(event, "cfs").particles();
 95
 96      // apply muon isolation
 97      Particles cand_mu;
 98      // pTcone around muon track
 99      for ( const Particle & mu :
100                apply<IdentifiedFinalState>(event,"muons").particlesByPt() ) {
101        double pTinCone = -mu.pT();
102        for ( const Particle & track : chg_tracks ) {
103          if ( deltaR(mu.momentum(),track.momentum()) < 0.2 )
104            pTinCone += track.pT();
105        }
106        if ( pTinCone < 1.8*GeV )
107          cand_mu.push_back(mu);
108      }
109
110      // Discard jets that overlap with electrons
111      Jets recon_jets;
112      for ( const Jet& jet : cand_jets ) {
113        bool away_from_e = true;
114        for ( const Particle & e : cand_e ) {
115          if ( deltaR(e.momentum(),jet.momentum()) <= 0.2 ) {
116            away_from_e = false;
117            break;
118          }
119        }
120        if ( away_from_e )
121         recon_jets.push_back( jet );
122      }
123
124      // Leptons far from jet
125      Particles recon_e;
126      for ( const Particle & e : cand_e ) {
127        bool e_near_jet = false;
128        for ( const Jet& jet : recon_jets ) {
129          if ( deltaR(e.momentum(),jet.momentum()) < 0.4 ) {
130            e_near_jet = true;
131            break;
132          }
133        }
134        // Electron isolation criterion
135        if ( ! e_near_jet ) {
136          double EtinCone = -e.Et();
137          for ( const Particle & track : chg_tracks) {
138            if ( deltaR(e.momentum(),track.momentum()) <= 0.2 )
139              EtinCone += track.Et();
140          }
141          if ( EtinCone/e.pT() <= 0.15 )
142            recon_e.push_back( e );
143        }
144      }
145
146      Particles recon_mu;
147      for ( const Particle & mu : cand_mu ) {
148         bool mu_near_jet = false;
149         for ( const Jet& jet : recon_jets ) {
150           if ( deltaR(mu.momentum(),jet.momentum()) < 0.4 ) {
151             mu_near_jet = true;
152             break;
153           }
154         }
155         if ( ! mu_near_jet )
156          recon_mu.push_back( mu );
157       }
158
159
160      // pTmiss
161      Particles vfs_particles
162        = apply<VisibleFinalState>(event, "vfs").particles();
163      FourMomentum pTmiss;
164      for ( const Particle & p : vfs_particles ) {
165        pTmiss -= p.momentum();
166      }
167      double eTmiss = pTmiss.pT();
168
169      // Exactly two leptons for each event
170      if ( recon_mu.size() + recon_e.size() != 2)
171        vetoEvent;
172
173      // Lepton pair mass
174      FourMomentum p_leptons;
175      for ( Particle e : recon_e ) {
176        p_leptons +=  e.momentum();
177      }
178      for ( Particle mu : recon_mu) {
179        p_leptons += mu.momentum();
180      }
181
182      if ( p_leptons.mass() <= 5.0 * GeV)
183        vetoEvent;
184
185      // ==================== FILL ====================
186
187
188      // electron, electron
189      if (recon_e.size() == 2 ) {
190
191        // SS ee
192        if ( recon_e[0].pid() * recon_e[1].pid() > 0 ) {
193          _hist_eTmiss_SS->fill(eTmiss);
194          if ( eTmiss > 100 ) {
195            MSG_DEBUG("Hits SS e+/-e+/-");
196            _count_SS_e_e->fill(0.5);
197          }
198        }
199
200        // OS ee
201        else if ( recon_e[0].pid() * recon_e[1].pid() < 0) {
202          _hist_eTmiss_OS->fill(eTmiss);
203          if ( eTmiss > 150 ) {
204            MSG_DEBUG("Hits OS e+e-");
205            _count_OS_e_e->fill(0.5);
206          }
207        }
208      }
209
210
211      // muon, electron
212      else if ( recon_e.size() == 1 ) {
213
214        // SS mu_e
215        if ( recon_e[0].pid() * recon_mu[0].pid() > 0 ) {
216          _hist_eTmiss_SS->fill(eTmiss);
217          if ( eTmiss > 100 ) {
218            MSG_DEBUG("Hits SS e+/-mu+/-");
219            _count_SS_e_mu->fill(0.5);
220          }
221        }
222
223        // OS mu_e
224        else if ( recon_e[0].pid() * recon_mu[0].pid() < 0) {
225          _hist_eTmiss_OS->fill(eTmiss);
226          if ( eTmiss > 150 ) {
227            MSG_DEBUG("Hits OS e+mu-");
228            _count_OS_e_mu->fill(0.5);
229          }
230        }
231      }
232
233
234      // muon, muon
235      else if ( recon_mu.size() == 2 ) {
236
237        // SS mu_mu
238        if ( recon_mu[0].pid() * recon_mu[1].pid() > 0 ) {
239          _hist_eTmiss_SS->fill(eTmiss);
240          if ( eTmiss > 100 ) {
241            MSG_DEBUG("Hits SS mu+/-mu+/-");
242            _count_SS_mu_mu->fill(0.5);
243          }
244        }
245
246        // OS mu_mu
247        else if ( recon_mu[0].pid() * recon_mu[1].pid() < 0) {
248          _hist_eTmiss_OS->fill(eTmiss);
249          if ( eTmiss > 150 ) {
250            MSG_DEBUG("Hits OS mu+mu-");
251            _count_OS_mu_mu->fill(0.5);
252          }
253        }
254      }
255
256
257    }
258
259    /// @}
260
261
262    void finalize() {
263      double norm = crossSection()/picobarn*35./sumOfWeights();
264      // event counts
265      scale(_count_OS_e_mu ,norm);
266      scale(_count_OS_e_e  ,norm);
267      scale(_count_OS_mu_mu,norm);
268      scale(_count_SS_e_mu ,norm);
269      scale(_count_SS_e_e  ,norm);
270      scale(_count_SS_mu_mu,norm);
271      scale(_hist_eTmiss_OS,10.*norm);
272      scale(_hist_eTmiss_SS,10.*norm);
273    }
274
275
276  private:
277
278    /// @name Histograms
279    /// @{
280    Histo1DPtr _count_OS_e_mu;
281    Histo1DPtr _count_OS_e_e;
282    Histo1DPtr _count_OS_mu_mu;
283    Histo1DPtr _count_SS_e_mu;
284    Histo1DPtr _count_SS_e_e;
285    Histo1DPtr _count_SS_mu_mu;
286    Histo1DPtr _hist_eTmiss_OS;
287    Histo1DPtr _hist_eTmiss_SS;
288    /// @}
289
290  };
291
292
293
294  RIVET_DECLARE_ALIASED_PLUGIN(ATLAS_2011_I894578, ATLAS_2011_S9019561);
295
296}