Rivet Analyses Reference

CDF_2004_S5839831

Transverse cone and `Swiss cheese' underlying event studies
Experiment: CDF (Tevatron Run 1)
Inspire ID: 647490
Status: VALIDATED
Authors:
  • Andy Buckley
References:Beams: p- p+
Beam energies: (315.0, 315.0); (900.0, 900.0) GeV
Run details:
  • QCD events at $\sqrt{s} = 630$ \& 1800 GeV. Several $p_\perp^\text{min}$ cutoffs are probably required to fill the profile histograms, e.g. 0 (min bias), 30, 90, 150 GeV at 1800 GeV, and 0 (min bias), 20, 90, 150 GeV at 630 GeV. Beam energy must be specified as analysis option "ENERGY" when rivet-merge'ing samples.

This analysis studies the underlying event via transverse cones of $R = 0.7$ at 90 degrees in \phi relative to the leading (highest E) jet, at $\sqrt{s} = 630$ and 1800 GeV. This is similar to the 2001 CDF UE analysis, except that cones, rather than the whole central \eta range are used. The transverse cones are categorised as `TransMIN' and `TransMAX' on an event-by-event basis, to give greater sensitivity to the UE component. `Swiss Cheese' distributions, where cones around the leading $n$ jets are excluded from the distributions, are also included for $n = 2, 3$. This analysis is useful for constraining the energy evolution of the underlying event, since it performs the same analyses at two distinct CoM energies. WARNING! The pT plots are normalised to raw number of events. The min bias data have not been reproduced by MC, and are not recommended for tuning. Beam energy must be specified (in GeV) as analysis option "ENERGY" when rivet-merging samples.

Source code: CDF_2004_S5839831.cc
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
// -*- C++ -*-
// "Acosta" underlying event analysis at CDF, inc. "Swiss Cheese"

#include "Rivet/Analysis.hh"
#include "Rivet/Jet.hh"
#include "Rivet/Projections/ChargedFinalState.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Projections/TriggerCDFRun0Run1.hh"

namespace Rivet {


  /// @brief CDF calo jet underlying event analysis at 630 and 1800 GeV
  ///
  /// CDF measurement of underlying event using calorimeter jet scales and
  /// alignment, particle flow activity in transverse cones, and the Swiss
  /// Cheese analysis method, where cones are excluded around the 2 and 3
  /// hardest jets.
  ///
  /// @author Andy Buckley
  class CDF_2004_S5839831 : public Analysis {
  public:

    RIVET_DEFAULT_ANALYSIS_CTOR(CDF_2004_S5839831);


    /// @name Analysis methods
    //@{

    void init() {
      // Set up projections
      declare(TriggerCDFRun0Run1(), "Trigger");
      const FinalState calofs(Cuts::abseta < 1.2);
      declare(calofs, "CaloFS");
      declare(FastJets(calofs, FastJets::CDFJETCLU, 0.7), "Jets");
      const ChargedFinalState trackfs(Cuts::abseta < 1.2 && Cuts::pT >= 0.4*GeV);
      declare(trackfs, "TrackFS");
      // Restrict tracks to |eta| < 0.7 for the min bias part.
      const ChargedFinalState mbfs(Cuts::abseta < 0.7 && Cuts::pT >= 0.4*GeV);
      declare(mbfs, "MBFS");
      // Restrict tracks to |eta| < 1 for the Swiss-Cheese part.
      const ChargedFinalState cheesefs(Cuts::abseta < 1.0 && Cuts::pT >= 0.4*GeV);
      declare(cheesefs, "CheeseFS");
      declare(FastJets(cheesefs, FastJets::CDFJETCLU, 0.7), "CheeseJets");

      // Book histograms
      if (isCompatibleWithSqrtS(1800)) {
        book(_pt90MaxAvg1800 ,1, 1, 1);
        book(_pt90MinAvg1800 ,1, 1, 2);
        book(_pt90Max1800 ,2, 1, 1);
        book(_pt90Min1800 ,2, 1, 2);
        book(_pt90Diff1800 ,2, 1, 3);
        book(_num90Max1800 ,4, 1, 1);
        book(_num90Min1800 ,4, 1, 2);
        book(_pTSum1800_2Jet ,7, 1, 1);
        book(_pTSum1800_3Jet ,7, 1, 2);

        book(_pt90Dbn1800Et40 ,3, 1, 1);
        book(_pt90Dbn1800Et80 ,3, 1, 2);
        book(_pt90Dbn1800Et120 ,3, 1, 3);
        book(_pt90Dbn1800Et160 ,3, 1, 4);
        book(_pt90Dbn1800Et200 ,3, 1, 5);
        book(_numTracksDbn1800MB ,5, 1, 1);
        book(_ptDbn1800MB ,6, 1, 1);
      } else if (isCompatibleWithSqrtS(630)) {
        book(_pt90Max630 ,8, 1, 1);
        book(_pt90Min630 ,8, 1, 2);
        book(_pt90Diff630 ,8, 1, 3);
        book(_pTSum630_2Jet ,9, 1, 1);
        book(_pTSum630_3Jet ,9, 1, 2);

        book(_numTracksDbn630MB ,10, 1, 1);
        book(_ptDbn630MB ,11, 1, 1);
      }
    }


    /// Do the analysis
    void analyze(const Event& event) {
      // Trigger
      const bool trigger = apply<TriggerCDFRun0Run1>(event, "Trigger").minBiasDecision();
      if (!trigger) vetoEvent;

      {
        MSG_DEBUG("Running max/min analysis");
        Jets jets = apply<JetAlg>(event, "Jets").jets(cmpMomByE);
        if (!jets.empty()) {
          // Leading jet must be in central |eta| < 0.5 region
          const Jet leadingjet = jets.front();
          const double etaLead = leadingjet.eta();
          // Get Et of the leading jet: used to bin histograms
          const double ETlead = leadingjet.Et();
          MSG_DEBUG("Leading Et = " << ETlead/GeV << " GeV");
          if (fabs(etaLead) > 0.5 && ETlead < 15*GeV) {
            MSG_DEBUG("Leading jet eta = " << etaLead
                     << " not in |eta| < 0.5 & pT > 15 GeV");
          } else {
            // Multiplicity & pT distributions for sqrt(s) = 630 GeV, 1800 GeV
            const Particles tracks = apply<FinalState>(event, "TrackFS").particles();
            const ConesInfo cones = _calcTransCones(leadingjet.momentum(), tracks);
            if (isCompatibleWithSqrtS(630)) {
              _pt90Max630->fill(ETlead/GeV, cones.ptMax/GeV);
              _pt90Min630->fill(ETlead/GeV, cones.ptMin/GeV);
              _pt90Diff630->fill(ETlead/GeV, cones.ptDiff/GeV);
            } else if (isCompatibleWithSqrtS(1800)) {
              _num90Max1800->fill(ETlead/GeV, cones.numMax);
              _num90Min1800->fill(ETlead/GeV, cones.numMin);
              _pt90Max1800->fill(ETlead/GeV, cones.ptMax/GeV);
              _pt90Min1800->fill(ETlead/GeV, cones.ptMin/GeV);
              _pt90Diff1800->fill(ETlead/GeV, cones.ptDiff/GeV);
              _pt90MaxAvg1800->fill(ETlead/GeV, cones.ptMax/GeV); // /numMax
              _pt90MinAvg1800->fill(ETlead/GeV, cones.ptMin/GeV); // /numMin
              //
              const double ptTransTotal = cones.ptMax + cones.ptMin;
              if (inRange(ETlead/GeV, 40., 80.)) {
                _pt90Dbn1800Et40->fill(ptTransTotal/GeV);
              } else if (inRange(ETlead/GeV, 80., 120.)) {
                _pt90Dbn1800Et80->fill(ptTransTotal/GeV);
              } else if (inRange(ETlead/GeV, 120., 160.)) {
                _pt90Dbn1800Et120->fill(ptTransTotal/GeV);
              } else if (inRange(ETlead/GeV, 160., 200.)) {
                _pt90Dbn1800Et160->fill(ptTransTotal/GeV);
              } else if (inRange(ETlead/GeV, 200., 270.)) {
                _pt90Dbn1800Et200->fill(ptTransTotal/GeV);
              }
            }

          }
        }
      }


      // Fill min bias total track multiplicity histos
      {
        MSG_DEBUG("Running min bias multiplicity analysis");
        const Particles mbtracks = apply<FinalState>(event, "MBFS").particles();
        if (isCompatibleWithSqrtS(1800)) {
          _numTracksDbn1800MB->fill(mbtracks.size());
        } else if (isCompatibleWithSqrtS(630)) {
          _numTracksDbn630MB->fill(mbtracks.size());
        }
        // Run over all charged tracks
        for (const Particle& t : mbtracks) {
          FourMomentum trackMom = t.momentum();
          const double pt = trackMom.pT();
          // Plot total pT distribution for min bias
          if (isCompatibleWithSqrtS(1800)) {
            _ptDbn1800MB->fill(pt/GeV);
          } else if (isCompatibleWithSqrtS(630)) {
            _ptDbn630MB->fill(pt/GeV);
          }
        }
      }



      // Construct "Swiss Cheese" pT distributions, with pT contributions from
      // tracks within R = 0.7 of the 1st, 2nd (and 3rd) jets being ignored. A
      // different set of charged tracks, with |eta| < 1.0, is used here, and all
      // the removed jets must have Et > 5 GeV.
      {
        MSG_DEBUG("Running Swiss Cheese analysis");
        const Particles cheesetracks = apply<FinalState>(event, "CheeseFS").particles();
        Jets cheesejets = apply<JetAlg>(event, "Jets").jets(cmpMomByE);
        if (cheesejets.empty()) {
          MSG_DEBUG("No 'cheese' jets found in event");
          return;
        }
        if (cheesejets.size() > 1 &&
            fabs(cheesejets[0].eta()) <= 0.5 &&
            cheesejets[0].Et()/GeV > 5.0 &&
            cheesejets[1].Et()/GeV > 5.0) {

          const double cheeseETlead = cheesejets[0].Et();

          const double eta1 = cheesejets[0].eta();
          const double phi1 = cheesejets[0].phi();
          const double eta2 = cheesejets[1].eta();
          const double phi2 = cheesejets[1].phi();

          double ptSumSub2(0), ptSumSub3(0);
          for (const Particle& t : cheesetracks) {
            FourMomentum trackMom = t.momentum();
            const double pt = trackMom.pT();

            // Subtracting 2 leading jets
            const double deltaR1 = deltaR(trackMom, eta1, phi1);
            const double deltaR2 = deltaR(trackMom, eta2, phi2);
            MSG_TRACE("Track vs jet(1): "
                     << "|(" << trackMom.eta() << ", " << trackMom.phi() << ") - "
                     << "|(" << eta1 << ", " << phi1 << ")| = " << deltaR1);
            MSG_TRACE("Track vs jet(2): "
                     << "|(" << trackMom.eta() << ", " << trackMom.phi() << ") - "
                     << "|(" << eta2 << ", " << phi2 << ")| = " << deltaR2);
            if (deltaR1 > 0.7 && deltaR2 > 0.7) {
              ptSumSub2 += pt;

              // Subtracting 3rd leading jet
              if (cheesejets.size() > 2 &&
                  cheesejets[2].Et()/GeV > 5.0) {
                const double eta3 = cheesejets[2].eta();
                const double phi3 = cheesejets[2].phi();
                const double deltaR3 = deltaR(trackMom, eta3, phi3);
                MSG_TRACE("Track vs jet(3): "
                         << "|(" << trackMom.eta() << ", " << trackMom.phi() << ") - "
                         << "|(" << eta3 << ", " << phi3 << ")| = " << deltaR3);
                if (deltaR3 > 0.7) {
                  ptSumSub3 += pt;
                }
              }
            }
          }

          // Swiss Cheese sub 2,3 jets distributions for sqrt(s) = 630 GeV, 1800 GeV
          if (isCompatibleWithSqrtS(630)) {
            if (!isZero(ptSumSub2)) _pTSum630_2Jet->fill(cheeseETlead/GeV, ptSumSub2/GeV);
            if (!isZero(ptSumSub3))_pTSum630_3Jet->fill(cheeseETlead/GeV, ptSumSub3/GeV);
          } else if (isCompatibleWithSqrtS(1800)) {
            if (!isZero(ptSumSub2))_pTSum1800_2Jet->fill(cheeseETlead/GeV, ptSumSub2/GeV);
            if (!isZero(ptSumSub3))_pTSum1800_3Jet->fill(cheeseETlead/GeV, ptSumSub3/GeV);
          }

        }
      }

    }


    void finalize() {
      /// @todo Take these normalisations from the data histo (it can't come from just the MC)

      if (isCompatibleWithSqrtS(1800)) {
        // Normalize to actual number of entries in pT dbn histos...
        normalize(_pt90Dbn1800Et40,  1656.75); // norm OK
        normalize(_pt90Dbn1800Et80,  4657.5); // norm OK
        normalize(_pt90Dbn1800Et120, 5395.5); // norm OK
        normalize(_pt90Dbn1800Et160, 7248.75); // norm OK
        normalize(_pt90Dbn1800Et200, 2442.0); // norm OK
      }

      // ...and for min bias distributions:
      if (isCompatibleWithSqrtS(1800)) {
        normalize(_numTracksDbn1800MB, 309718.25); // norm OK
        normalize(_ptDbn1800MB, 33600.0); // norm OK
      } else if (isCompatibleWithSqrtS(630)) {
        normalize(_numTracksDbn630MB, 1101024.0); // norm OK
        normalize(_ptDbn630MB, 105088.0); // norm OK
      }
    }

    //@}


  private:


    /// @name Cone machinery
    /// @{

    /// @cond CONEUE_DETAIL

    struct ConesInfo {
      ConesInfo() : numMax(0), numMin(0), ptMax(0), ptMin(0), ptDiff(0) {}
      unsigned int numMax, numMin;
      double ptMax, ptMin, ptDiff;
    };

    /// @endcond


    ConesInfo _calcTransCones(const double etaLead, const double phiLead,
                              const Particles& tracks) {
      const double phiTransPlus = mapAngle0To2Pi(phiLead + PI/2.0);
      const double phiTransMinus = mapAngle0To2Pi(phiLead - PI/2.0);
      MSG_DEBUG("phi_lead = " << phiLead
               << " -> trans = (" << phiTransPlus
               << ", " << phiTransMinus << ")");

      unsigned int numPlus(0), numMinus(0);
      double ptPlus(0), ptMinus(0);
      // Run over all charged tracks
      for (const Particle& t : tracks) {
        FourMomentum trackMom = t.momentum();
        const double pt = trackMom.pT();

        // Find if track mom is in either transverse cone
        if (deltaR(trackMom, etaLead, phiTransPlus) < 0.7) {
          ptPlus += pt;
          numPlus += 1;
        } else if (deltaR(trackMom, etaLead, phiTransMinus) < 0.7) {
          ptMinus += pt;
          numMinus += 1;
        }
      }

      ConesInfo rtn;
      // Assign N_{min,max} from N_{plus,minus}
      rtn.numMax = (ptPlus >= ptMinus) ? numPlus : numMinus;
      rtn.numMin = (ptPlus >= ptMinus) ? numMinus : numPlus;
      // Assign pT_{min,max} from pT_{plus,minus}
      rtn.ptMax = (ptPlus >= ptMinus) ? ptPlus : ptMinus;
      rtn.ptMin = (ptPlus >= ptMinus) ? ptMinus : ptPlus;
      rtn.ptDiff = fabs(rtn.ptMax - rtn.ptMin);

      MSG_DEBUG("Min cone has " << rtn.numMin << " tracks -> "
               << "pT_min = " << rtn.ptMin/GeV << " GeV");
      MSG_DEBUG("Max cone has " << rtn.numMax << " tracks -> "
               << "pT_max = " << rtn.ptMax/GeV << " GeV");

      return rtn;
    }


    ConesInfo _calcTransCones(const FourMomentum& leadvec,
                              const Particles& tracks) {
      const double etaLead = leadvec.eta();
      const double phiLead = leadvec.phi();
      return _calcTransCones(etaLead, phiLead, tracks);
    }

    /// @}


    /// @name Histogram collections
    /// @{

    /// Profile histograms, binned in the \f$ E_T \f$ of the leading jet, for
    /// the average \f$ p_T \f$ in the toward, transverse and away regions at
    /// \f$ \sqrt{s} = 1800 \text{GeV} \f$.
    /// Corresponds to Table 1, and HepData table 1.
    Profile1DPtr _pt90MaxAvg1800, _pt90MinAvg1800;

    /// Profile histograms, binned in the \f$ E_T \f$ of the leading jet, for
    /// the \f$ p_T \f$ sum in the toward, transverse and away regions at
    /// \f$ \sqrt{s} = 1800 \text{GeV} \f$.
    /// Corresponds to figure 2/3, and HepData table 2.
    Profile1DPtr _pt90Max1800, _pt90Min1800, _pt90Diff1800;

    /// Profile histograms, binned in the \f$ E_T \f$ of the leading jet, for
    /// the \f$ p_T \f$ sum in the toward, transverse and away regions at
    /// at \f$ \sqrt{s} = 630 \text{GeV} \f$.
    /// Corresponds to figure 8, and HepData table 8.
    Profile1DPtr _pt90Max630, _pt90Min630, _pt90Diff630;

    /// Profile histograms, binned in the \f$ E_T \f$ of the leading jet, for
    /// the cone track multiplicity at \f$ \sqrt{s} = 1800 \text{GeV} \f$.
    /// Corresponds to figure 5, and HepData table 4.
    Profile1DPtr _num90Max1800, _num90Min1800;

    /// Profile histograms, binned in the \f$ E_T \f$ of the leading jet, for
    /// the \f$ p_T \f$ sum at \f$ \sqrt{s} = 1800 \text{GeV} \f$.
    /// Corresponds to figure 7, and HepData table 7.
    Profile1DPtr _pTSum1800_2Jet, _pTSum1800_3Jet;

    /// Profile histograms, binned in the \f$ E_T \f$ of the leading jet, for
    /// the \f$ p_T \f$ sum at \f$ \sqrt{s} = 630 \text{GeV} \f$.
    /// Corresponds to figure 9, and HepData table 9.
    Profile1DPtr _pTSum630_2Jet, _pTSum630_3Jet;

    /// Histogram of \f$ p_{T\text{sum}} \f$ distribution for 5 different
    /// \f$ E_{T1} \f$ bins.
    /// Corresponds to figure 4, and HepData table 3.
    Histo1DPtr _pt90Dbn1800Et40, _pt90Dbn1800Et80, _pt90Dbn1800Et120,
      _pt90Dbn1800Et160, _pt90Dbn1800Et200;

    /// Histograms of track multiplicity and \f$ p_T \f$ distributions for
    /// minimum bias events.
    /// Figure 6, and HepData tables 5 & 6.
    /// Figure 10, and HepData tables 10 & 11.
    Histo1DPtr _numTracksDbn1800MB, _ptDbn1800MB;
    Histo1DPtr _numTracksDbn630MB, _ptDbn630MB;

    /// @}

  };



  RIVET_DECLARE_ALIASED_PLUGIN(CDF_2004_S5839831, CDF_2004_I647490);

}