Rivet Analyses Reference

ATLAS_2012_I1084540

Rapidity gap cross sections measured with the ATLAS detector in $pp$ collisions at $\sqrt{s} = 7$ TeV.
Experiment: ATLAS (LHC 7TeV)
Inspire ID: 1084540
Status: VALIDATED
Authors:
  • Oldrich Kepka
  • Tim Martin
  • Paul Newman
  • Pavel Ruzicka
References:Beams: p+ p+
Beam energies: (3500.0, 3500.0) GeV
Run details:
  • Minimum bias inelastic pp collision at 7 TeV including diffractive component and overall cross section.

Pseudorapidity gap distributions in proton-proton collisions at $\sqrt{s} = 7$ TeV are studied using a minimum bias data sample with an integrated luminosity of 7.1 inverse microbarns. Cross sections are measured differentially in terms of $\Delta \eta_F$, the larger of the pseudorapidity regions extending to the limits of the ATLAS sensitivity, at $\eta = \pm 4.9$, in which no final state particles are produced above a transverse momentum threshold $p_\perp$ cut. The measurements span the region $0 < \Delta \eta_F < 8$ for $200 < p_T\text{ cut} < 800 \text{MeV}$. At small $\Delta \eta_F$, the data test the reliability of hadronisation models in describing rapidity and transverse momentum fluctuations in final state particle production. The measurements at larger gap sizes are dominated by contributions from the single diffractive dissociation process ($pp \to Xp$), enhanced by double dissociation ($pp \to XY$) where the invariant mass of the lighter of the two dissociation systems satisfies $M_Y \lesssim 7 \text{GeV}$. The resulting cross section is $\mathrm{d} \sigma / \mathrm{d} \Delta \eta_F \sim 1$ mb for $\Delta \eta_F \gtrsim 3$. The large rapidity gap data are used to constrain the value of the pomeron intercept appropriate to triple Regge models of soft diffraction. The cross section integrated over all gap sizes is compared with other LHC inelastic cross section measurements.

Source code: ATLAS_2012_I1084540.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

// -*- C++ -*-
/**
 * @name ATLAS Diffractive Gaps Rivet Analysis
 * @author Tim Martin, tim.martin@cern.ch
 * @version 1.0
 * @date 16/01/2012
 * @see http://arxiv.org/abs/1201.2808
 * @note pp, sqrt(s) = 7 TeV
 * @note Rapidity gap finding algorithm designed to compliment
 * the ATLAS detector acceptance. Forward rapidity gap sizes
 * are calculated for each event, considering all stable
 * particles above pT cut values 200, 400, 600 and 800 MeV in
 * turn. A forward rapidity gap is defined to be the largest
 * continuous region stretching inward from either edge of the
 * detector at eta = |4.9| which contains zero particles above
 * pT Cut. Soft diffractive topologies are isolated at large
 * gap sizes.
 *
 */
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"

namespace Rivet {


  class ATLAS_2012_I1084540 : public Analysis {
  public:

    ATLAS_2012_I1084540() : Analysis("ATLAS_2012_I1084540") {}


    /// @name Analysis methods
    //@{
    /// Book histograms and initialise projections before the run
    void init() {
      //All final states. Rapidity range = ATLAS calorimetry. Lowest pT cut = 200 MeV.
      const FinalState cnfs2((Cuts::etaIn(-_etaMax, _etaMax) && Cuts::pT >=  0.2 * GeV));
      const FinalState cnfs4((Cuts::etaIn(-_etaMax, _etaMax) && Cuts::pT >=  0.4 * GeV));
      const FinalState cnfs6((Cuts::etaIn(-_etaMax, _etaMax) && Cuts::pT >=  0.6 * GeV));
      const FinalState cnfs8((Cuts::etaIn(-_etaMax, _etaMax) && Cuts::pT >=  0.8 * GeV));
      declare(cnfs2, "CNFS2");
      declare(cnfs4, "CNFS4");
      declare(cnfs6, "CNFS6");
      declare(cnfs8, "CNFS8");

      _etaBinSize = (2. * _etaMax)/(double)_etaBins;

      //Book histogram
      book(_h_DeltaEtaF_200 ,1, 1, 1);
      book(_h_DeltaEtaF_400 ,2, 1, 1);
      book(_h_DeltaEtaF_600 ,3, 1, 1);
      book(_h_DeltaEtaF_800 ,4, 1, 1);
    }

  private:
    void fillMap(const FinalState& fs, bool* energyMap, double pTcut) {
      // Fill true/false array by iterating over particles and compare their
      // pT with pTcut
      for (const Particle& p : fs.particles(cmpMomByEta)) {
        int checkBin = -1;
        double checkEta = -_etaMax;
        while (1) {
          checkEta += _etaBinSize;
          ++checkBin;
          if (p.eta() < checkEta) {
            energyMap[checkBin] = (p.pT() > pTcut * GeV);
            break;
          }
        }
      }
    }

  public:
    /// Perform the per-event analysis
    void analyze(const Event& event) {
      static unsigned int event_count = 0;
      ++event_count;
      const FinalState& fs2 = apply<FinalState>(event, "CNFS2");
      const FinalState& fs4 = apply<FinalState>(event, "CNFS4");
      const FinalState& fs6 = apply<FinalState>(event, "CNFS6");
      const FinalState& fs8 = apply<FinalState>(event, "CNFS8");

      // Set up Yes/No arrays for energy in each eta bin at each pT cut
      bool energyMap_200[_etaBins];
      bool energyMap_400[_etaBins];
      bool energyMap_600[_etaBins];
      bool energyMap_800[_etaBins];
      for (int i = 0; i < _etaBins; ++i) {
        energyMap_200[i] = false;
        energyMap_400[i] = false;
        energyMap_600[i] = false;
        energyMap_800[i] = false;
      }

      // Veto bins based on final state particles > Cut (Where Cut = 200 - 800 MeV pT)
      fillMap(fs2, energyMap_200, 0.2);
      fillMap(fs4, energyMap_400, 0.4);
      fillMap(fs6, energyMap_600, 0.6);
      fillMap(fs8, energyMap_800, 0.8);

      // Apply gap finding algorithm
      // Detector layout follows...
      // -Eta [Proton  ---  DetectorCSide  ---  DetectorBarrel  ---  DetectorASide  ---  Proton] +Eta
      bool gapDirectionAt200 = false; //False is gap on C size, True is gap on A side.
      double largestEdgeGap_200 = 0.;
      double largestEdgeGap_400 = 0.;
      double largestEdgeGap_600 = 0.;
      double largestEdgeGap_800 = 0.;

      for (int E = 200; E <= 800; E += 200) {
        double EdgeGapSizeA = -1, EdgeGapSizeC = -1;
        bool* energyMap = 0;
        switch (E) {
          case 200: energyMap = energyMap_200; break;
          case 400: energyMap = energyMap_400; break;
          case 600: energyMap = energyMap_600; break;
          case 800: energyMap = energyMap_800; break;
        }

        // Look left to right
        for (int a = 0; a < _etaBins; ++a) {
          if (energyMap[a] == true) {
            EdgeGapSizeA = (_etaBinSize * a);
            break;
          }
        }

        // And look right to left
        for (int c = _etaBins-1; c >= 0; --c) {
          if (energyMap[c] == true) {
            EdgeGapSizeC = (2 * _etaMax) - (_etaBinSize * (c+1));
            if (fuzzyEquals(EdgeGapSizeC, 4.47035e-08)) EdgeGapSizeC = 0.0;
            break;
          }
        }
        // Put your hands on your hips

        // Find the largest gap
        double largestEdgeGap = 0.;
        if (E == 200) {
          // If the 200 MeV pass, take the biggest of the two gaps. Make note of which side for higher pT cuts.
          largestEdgeGap = std::max(EdgeGapSizeA,EdgeGapSizeC);
          gapDirectionAt200 = (EdgeGapSizeA > EdgeGapSizeC);
        } else {
          // Use the direction from 200 MeV pass, most accurate measure of which side gap is on.
          if (gapDirectionAt200) {
            largestEdgeGap = EdgeGapSizeA;
          }
          else largestEdgeGap = EdgeGapSizeC;
        }

        // Check case of empty detector
        if (largestEdgeGap < 0.0) largestEdgeGap = 2.0 * _etaMax;

        // Fill bin centre
        switch (E) {
          case 200: _h_DeltaEtaF_200->fill(largestEdgeGap + _etaBinSize/2.); break;
          case 400: _h_DeltaEtaF_400->fill(largestEdgeGap + _etaBinSize/2.); break;
          case 600: _h_DeltaEtaF_600->fill(largestEdgeGap + _etaBinSize/2.); break;
          case 800: _h_DeltaEtaF_800->fill(largestEdgeGap + _etaBinSize/2.); break;
        }

        if (E == 200) largestEdgeGap_200 = largestEdgeGap;
        if (E == 400) largestEdgeGap_400 = largestEdgeGap;
        if (E == 600) largestEdgeGap_600 = largestEdgeGap;
        if (E == 800) largestEdgeGap_800 = largestEdgeGap;
      }

      // Algorithm result every 1000 events
      if (event_count % 1000 == 0) {
        for (int E = 200; E <= 800; E += 200) {
          bool* energyMap = 0;
          double largestEdgeGap = 0;
          switch (E) {
            case 200: energyMap = energyMap_200; largestEdgeGap = largestEdgeGap_200; break;
            case 400: energyMap = energyMap_400; largestEdgeGap = largestEdgeGap_400; break;
            case 600: energyMap = energyMap_600; largestEdgeGap = largestEdgeGap_600; break;
            case 800: energyMap = energyMap_800; largestEdgeGap = largestEdgeGap_800; break;
          }
          MSG_DEBUG("Largest Forward Gap at pT Cut " << E << " MeV=" << largestEdgeGap
            << " eta, NFinalState pT > 200 in ATLAS acceptance:" << fs2.particles().size());
          std::string hitPattern = "Detector HitPattern=-4.9[";
          for (int a = 0; a < _etaBins; ++a) {
            if (energyMap[a] == true) hitPattern += "X";
            else hitPattern += "_";
          }
          hitPattern += "]4.9";
          MSG_DEBUG(hitPattern);
          std::string gapArrow = "                         ";
            if (!gapDirectionAt200) {
            int drawSpaces = (int)(_etaBins - (largestEdgeGap/_etaBinSize) + 0.5);
            for (int i = 0; i < drawSpaces; ++i) gapArrow += " ";
          }
          int drawArrows = (int)((largestEdgeGap/_etaBinSize) + 0.5);
          for (int i = 0; i < drawArrows; ++i) gapArrow += "^";
          MSG_DEBUG(gapArrow);
        }
      }
    }

    /// Normalise histograms after the run, Scale to cross section
    void finalize() {
      MSG_DEBUG("Cross Section=" << crossSection() / millibarn << "mb, SumOfWeights=" << sumOfWeights());
      scale(_h_DeltaEtaF_200, (crossSection() / millibarn)/sumOfWeights());
      scale(_h_DeltaEtaF_400, (crossSection() / millibarn)/sumOfWeights());
      scale(_h_DeltaEtaF_600, (crossSection() / millibarn)/sumOfWeights());
      scale(_h_DeltaEtaF_800, (crossSection() / millibarn)/sumOfWeights());
    }
    //@}

  private:
    /// @name Histograms
    //@{
    Histo1DPtr _h_DeltaEtaF_200;
    Histo1DPtr _h_DeltaEtaF_400;
    Histo1DPtr _h_DeltaEtaF_600;
    Histo1DPtr _h_DeltaEtaF_800;
    //@}
    /// @name Private variables
    //@{
    static constexpr int _etaBins = 49;
    static constexpr double _etaMax = 4.9;
    double _etaBinSize;
    //@}
  };

  // The hook for the plugin system
  RIVET_DECLARE_PLUGIN(ATLAS_2012_I1084540);

}

Edit this page on GitLab