rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

H1_1996_I424463

Transverse momentum spectra of charged particles in DIS (H1 1996)
Experiment: H1 (HERA)
Inspire ID: 424463
Status: VALIDATED
Authors:
  • Suraj Kumar Singh
  • Hannes Jung
  • Andrii Verbytskyi
References: Beams: e+ p+, p+ e+
Beam energies: (27.5, 820.0); (820.0, 27.5) GeV
    No run details listed

Transverse momentum spectra of charged particles produced in deep inelastic scattering are measured as a function of the kinematic variables x and Q2 using the H1 detector at the epcollider HERA. The data are compared to different patton emission models, either with or without ordering of the emissions in transverse momentum. The data provide evidence for a relatively large amount of patton radiation between the current and the remnant systems.

Source code: H1_1996_I424463.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Tools/ParticleIdUtils.hh"
  4#include "Rivet/Projections/FinalState.hh"
  5#include "Rivet/Projections/DISKinematics.hh"
  6#include "Rivet/Projections/ChargedFinalState.hh"
  7#include "Rivet/Projections/DISKinematics.hh"
  8
  9namespace Rivet {
 10
 11
 12  /// @brief Transverse momentum spectra of charged particles in DIS (H1 1996)
 13  class H1_1996_I424463 : public Analysis {
 14  public:
 15
 16    /// Constructor
 17    DEFAULT_RIVET_ANALYSIS_CTOR(H1_1996_I424463);
 18
 19
 20    /// @name Analysis methods
 21    ///@{
 22
 23    /// Book histograms and initialise projections before the run
 24    void init() {
 25
 26      // Book histograms
 27      declare(DISLepton(), "Lepton");
 28      declare(DISKinematics(), "Kinematics");
 29      declare(ChargedFinalState(), "CFS");
 30      declare(FinalState(), "FS");
 31
 32
 33      book(_NevAll, "TMP/Nev_all");
 34      _h_dndpt_high_eta_bin.resize(10);
 35      _h_dndpt_low_eta_bin.resize(10);
 36      _hdndeta_pt1_bin.resize(10);
 37      _hdndeta_bin.resize(10);
 38      _hdndptmax_low_eta_bin.resize(8);
 39      int ixx = 0 ;
 40      for (size_t ix = 0; ix < 10; ++ix) {
 41        book(_Nevt_after_cuts[ix], "TMP/Nevt_after_cuts" + to_string(ix));
 42        book(_h_dndpt_high_eta_bin[ix], ix+1, 1, 1);
 43        book(_h_dndpt_low_eta_bin[ix], ix+11, 1, 1);
 44        if (ix != 6 && ix != 9) {
 45          book(_hdndptmax_low_eta_bin[ixx], ixx+21, 1, 1);
 46          ixx=ixx+1;
 47        }
 48        book(_hdndeta_pt1_bin[ix], ix+29, 1, 1);
 49        book(_hdndeta_bin[ix],  ix+39, 1, 1);
 50      }
 51
 52
 53    }
 54
 55
 56    /// Perform the per-event analysis
 57    void analyze(const Event& event) {
 58      // const ChargedFinalState& cfs = applyProjection<ChargedFinalState>(event, "CFS");
 59      const FinalState& fs = apply<FinalState>(event, "FS");
 60      const DISKinematics& dk = apply<DISKinematics>(event, "Kinematics");
 61      const DISLepton& dl = apply<DISLepton>(event,"Lepton");
 62
 63      // Get the DIS kinematics
 64      double x  = dk.x();
 65      double y = dk.y();
 66      double Q2 = dk.Q2()/GeV;
 67      double W2 = dk.W2()/GeV;
 68
 69      // Momentum of the scattered lepton
 70      FourMomentum leptonMom = dl.out().momentum();
 71      double enel = leptonMom.E()/GeV;
 72      double thel = 180.-leptonMom.angle(dl.in().momentum())/degree;
 73
 74      _NevAll -> fill() ;
 75
 76      const bool cut = y > 0.05 && Q2 > 5. && Q2 < 100.&& enel > 12. && W2 > 4400. && thel > 157. && thel < 173.;
 77      if (!cut) vetoEvent;
 78
 79      int ibin[10] ;
 80      for ( int i=0; i< 10; i++) {
 81      ibin[i] = 0; }
 82
 83      if(5.<Q2&&Q2<50.&&x>0.0001&&x<0.001)    ibin[0]=1;
 84      if(5.<Q2&&Q2<10.&&x>0.0001&&x<0.0002)   ibin[1]=1;
 85      if(6.<Q2&&Q2<10.&&x>0.0002&&x<0.0005)   ibin[2]=1;
 86      if(10.<Q2&&Q2<20.&&x>0.0002&&x<0.0005)  ibin[3]=1;
 87      if(10.<Q2&&Q2<20.&&x>0.0005&&x<0.0008)  ibin[4]=1;
 88      if(10.<Q2&&Q2<20.&&x>0.0008&&x<0.0015)  ibin[5]=1;
 89      if(10.<Q2&&Q2<20.&&x>0.0015&&x<0.004)   ibin[6]=1;
 90      if(20.<Q2&&Q2<50.&&x>0.0005&&x<0.0014)  ibin[7]=1;
 91      if(20.<Q2&&Q2<50.&&x>0.0014&&x<0.003)   ibin[8]=1;
 92      if(20.<Q2&&Q2<50.&&x>0.003&&x<0.01)     ibin[9]=1;
 93
 94      for ( int i=0; i< 10; i++) { if(ibin[i]==1) _Nevt_after_cuts[i] ->fill(); }
 95
 96      // Extract the particles other than the lepton
 97      Particles particles;
 98      particles.reserve(fs.particles().size());
 99      ConstGenParticlePtr dislepGP = dl.out().genParticle();
100      for (const Particle& p : fs.particles()) {
101        ConstGenParticlePtr loopGP = p.genParticle();
102        if (loopGP == dislepGP) continue;
103        particles.push_back(p);
104      }
105
106      // Boost to hadronic CM
107      const LorentzTransform hcmboost = dk.boostHCM();
108
109      int mult = 0 ;
110      // Loop over the particles
111      // long ncharged(0);
112      double ptmax_high[10], ptmax_low[10] ;
113      for ( int i=0; i< 10; i++) {
114       ptmax_high[i] = 0.; ptmax_low[i] = 0.; }
115      double EtSum = 0;
116      for (size_t ip1 = 0; ip1 < particles.size(); ++ip1) {
117         const Particle& p = particles[ip1];
118
119         double eta = p.momentum().pseudorapidity();
120         // Boost to hcm
121         const FourMomentum hcmMom = hcmboost.transform(p.momentum());
122
123         // apply safety cuts
124         if (eta > -5  && eta < 10.) {
125            mult = mult + 1;
126            double pThcm =hcmMom.pT() ;
127            double etahcm = hcmMom.pseudorapidity();
128            // cout << " charge " << PID::charge(p.pid()) << endl;
129            if (etahcm > 0. && etahcm < 2.0) {
130              EtSum = EtSum + hcmMom.Et();
131            }
132            if (PID::charge(p.pid()) != 0) {
133              if (etahcm > 0.5 && etahcm < 1.5) {
134               for ( int i=0; i< 10; i++) {
135                 if (ibin[i]==1) {
136                   _h_dndpt_low_eta_bin[i] ->fill(pThcm);
137                   if (pThcm > ptmax_low[i] ) ptmax_low[i] = pThcm;
138                 }
139               }
140
141              }
142              if (etahcm > 1.5 && etahcm < 2.5){
143               for ( int i=0; i< 10; i++) {
144                 if (ibin[i]==1) {
145                   _h_dndpt_high_eta_bin[i] ->fill(pThcm);
146                   if (pThcm > ptmax_high[i] ) ptmax_high[i] = pThcm;
147                 }
148               }
149              }
150              for ( int i=0; i< 10; ++i) {
151                if(ibin[i]==1) _hdndeta_bin[i] ->fill(etahcm);
152                if(ibin[i]==1 && pThcm > 1.) _hdndeta_pt1_bin[i] ->fill(etahcm);
153              }
154           }
155         }  // end of loop over the particles
156      }
157      int ii=0;
158      for ( int i=0; i< 10; ++i) {
159         if (i != 6 && i != 9 ) {
160             if( ibin[i]==1 && EtSum > 6. ) {
161              _hdndptmax_low_eta_bin[ii] ->fill(ptmax_low[i]);
162             }
163            ii=ii+1;
164         }
165      }
166
167    }
168    /// Normalise histograms etc., after the run
169    void finalize() {
170
171
172      int ii =0;
173      for (int i=0; i < 10; ++i) {
174        if (_Nevt_after_cuts[i]->val() != 0) {
175          scale(_h_dndpt_high_eta_bin[i], 1./ *_Nevt_after_cuts[i]);
176          scale(_h_dndpt_low_eta_bin[i],  1./ *_Nevt_after_cuts[i]);
177          scale(_hdndeta_bin[i],          1./ *_Nevt_after_cuts[i]);
178          scale(_hdndeta_pt1_bin[i],      1./ *_Nevt_after_cuts[i]); 
179        }
180        if ( i !=6 && i!=9) {
181          if (_Nevt_after_cuts[i]->val()  != 0) normalize(_hdndptmax_low_eta_bin[ii]); ii=ii+1; 
182        }
183      }
184    }
185
186    ///@}
187
188
189    /// @name Histograms
190    ///@{
191    map<string, Histo1DPtr> _h;
192    map<string, Profile1DPtr> _p;
193    map<string, CounterPtr> _c;
194    array<CounterPtr,10> _Nevt_after_cuts;
195    vector<Histo1DPtr> _h_dndpt_low_eta_bin, _h_dndpt_high_eta_bin;
196    vector<Histo1DPtr> _hdndeta_bin, _hdndeta_pt1_bin, _hdndptmax_low_eta_bin;
197    CounterPtr _NevAll;
198    ///@}
199
200
201  };
202
203
204  RIVET_DECLARE_ALIASED_PLUGIN(H1_1996_I424463, H1_1997_I424463);
205
206}