rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

H1_2007_I746380

Tests of QCD Factorisation in the diffractive production of dijets in deep-inelastic scattering and photoproduction at HERA
Experiment: H1 (HERA)
Inspire ID: 746380
Status: UNVALIDATED
Authors:
  • Christine O. Rasmussen
  • Ilkka Helenius
References: Beams: p+ e+
Beam energies: (820.0, 27.5) GeV
Run details:
  • 820 GeV protons colliding with 27.5 GeV positrons; Diffractive photoproduction of dijets; Jet $E_T > 4$ GeV;

H1 diffractive jets from proton--positron collisions at beam energies of 820 GeV and 27.5 GeV. Measurements are presented of differential dijet cross sections in diffractive photoproduction (Q2 < 0.01 GeV2) and deep-inelastic scattering processes (DIS, 4 < Q2 < 80 GeV2). The event topology is given by ep → eXY , in which the system X, containing at least two jets, is separated from a leading low-mass proton remnant system Y by a large rapidity gap.

Source code: H1_2007_I746380.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/Beam.hh"
  4#include "Rivet/Projections/FinalState.hh"
  5#include "Rivet/Projections/DISKinematics.hh"
  6#include "Rivet/Projections/DISFinalState.hh"
  7#include "Rivet/Projections/FastJets.hh"
  8
  9namespace Rivet {
 10
 11  namespace H1_2007_I746380_PROJECTIONS {
 12
 13    /// Projection to find the largest gaps and the masses of the two
 14    /// systems separated by the gap. Based on the HZTools gap-finding
 15    /// method (hzhadgap.F). Note that gaps are found in the HCM frame.
 16    ///
 17    /// @author Christine O. Rasmussen.
 18    class RapidityGap : public Projection {
 19
 20    public:
 21
 22      /// Type of DIS boost to apply
 23      enum Frame { HCM, LAB, XCM };
 24
 25      RapidityGap() {
 26        setName("RapidityGap");
 27        declare(DISKinematics(), "DISKIN");
 28        declare(DISFinalState(DISFinalState::BoostFrame::HCM), "DISFS");
 29      }
 30
 31      DEFAULT_RIVET_PROJ_CLONE(RapidityGap);
 32
 33      double M2X() const { return _M2X; }
 34      double M2Y() const { return _M2Y; }
 35      double t() const { return _t; }
 36      double gap() const { return _gap; }
 37      double gapUpp() const { return _gapUpp; }
 38      double gapLow() const { return _gapLow; }
 39      double EpPzX(Frame f) const {
 40        if (f == LAB) return _ePpzX_LAB;
 41        else if (f == XCM) return _ePpzX_XCM;
 42        else return _ePpzX_HCM;
 43      }
 44      double EmPzX(Frame f) const {
 45        if (f == LAB) return _eMpzX_LAB;
 46        else if (f == XCM) return _eMpzX_XCM;
 47        else return _eMpzX_HCM;
 48      }
 49      FourMomentum pX(Frame f) const {
 50        if (f == LAB) return _momX_LAB;
 51        else if (f == XCM) return _momX_XCM;
 52        else return _momX_HCM;
 53      }
 54      FourMomentum pY(Frame f) const {
 55        if (f == LAB) return _momY_LAB;
 56        else if (f == XCM) return _momY_XCM;
 57        else return _momY_HCM;
 58      }
 59      const Particles& systemX(Frame f) const {
 60        if (f == LAB) return _pX_LAB;
 61        else if (f == XCM) return _pX_XCM;
 62        else return _pX_HCM;
 63      }
 64      const Particles& systemY(Frame f) const {
 65        if (f == LAB) return _pY_LAB;
 66        else if (f == XCM) return _pY_XCM;
 67        else return _pY_HCM;
 68      }
 69
 70    protected:
 71
 72      virtual CmpState compare(const Projection& p) const {
 73        const RapidityGap& other = pcast<RapidityGap>(p);
 74        return mkNamedPCmp(other, "DISKIN") || mkNamedPCmp(other, "DISFS");
 75      }
 76
 77      virtual void project(const Event& e){
 78        const DISKinematics& dk = apply<DISKinematics>(e, "DISKIN");
 79        const Particles& p      = apply<DISFinalState>(e, "DISFS").particles(cmpMomByEta);
 80        findgap(p, dk);
 81      }
 82
 83      void clearAll(){
 84        _M2X = _M2Y = _t = _gap = 0.;
 85        _gapUpp = _gapLow = -8.;
 86        _ePpzX_HCM = _eMpzX_HCM =_ePpzX_LAB = _eMpzX_LAB = _ePpzX_XCM = _eMpzX_XCM = 0.;
 87        _momX_HCM.setPE(0., 0., 0., 0.);
 88        _momY_HCM.setPE(0., 0., 0., 0.);
 89        _momX_XCM.setPE(0., 0., 0., 0.);
 90        _momY_XCM.setPE(0., 0., 0., 0.);
 91        _momX_LAB.setPE(0., 0., 0., 0.);
 92        _momY_LAB.setPE(0., 0., 0., 0.);
 93        _pX_HCM.clear();
 94        _pY_HCM.clear();
 95        _pX_XCM.clear();
 96        _pY_XCM.clear();
 97        _pX_LAB.clear();
 98        _pY_LAB.clear();
 99      }
100
101      void findgap(const Particles& particles, const DISKinematics& diskin){
102
103        clearAll();
104
105        // Begin by finding largest gap and gapedges between all final
106        // state particles in HCM frame.
107        int nP  = particles.size();
108        int dir = diskin.orientation();
109        for (int i = 0; i < nP-1; ++i){
110          double tmpGap = abs(particles[i+1].eta() - particles[i].eta());
111          if (tmpGap > _gap) {
112            _gap    = tmpGap;
113            _gapLow = (dir > 0) ? particles[i].eta() : dir * particles[i+1].eta();
114            _gapUpp = (dir > 0) ? particles[i+1].eta() : dir * particles[i].eta();
115          }
116        }
117
118        // Define the two systems X and Y.
119        Particles tmp_pX, tmp_pY;
120        for (const Particle& ip : particles) {
121          if (dir * ip.eta() > _gapLow) tmp_pX.push_back(ip);
122          else tmp_pY.push_back(ip);
123        }
124
125        Particles pX, pY;
126        pX = (dir < 0) ? tmp_pY : tmp_pX;
127        pY = (dir < 0) ? tmp_pX : tmp_pY;
128
129        // Find variables related to HCM frame.
130        // Note that HCM has photon along +z, as opposed to
131        // H1 where proton is along +z. This results in a sign change
132        // as compared to H1 papers!
133
134        // X - side
135        FourMomentum momX;
136        for (const Particle& jp : pX) {
137          momX  += jp.momentum();
138          _ePpzX_HCM += jp.E() - jp.pz(); // Sign + => -
139          _eMpzX_HCM += jp.E() + jp.pz(); // Sign - => +
140        }
141        _momX_HCM = momX;
142        _pX_HCM   = pX;
143        _M2X      = _momX_HCM.mass2();
144
145        // Y - side
146        FourMomentum momY;
147        for (const Particle& kp : pY) momY += kp.momentum();
148        _momY_HCM = momY;
149        _pY_HCM   = pY;
150        _M2Y      = _momY_HCM.mass2();
151
152        // Find variables related to LAB frame
153        const LorentzTransform hcmboost   = diskin.boostHCM();
154        const LorentzTransform hcminverse = hcmboost.inverse();
155        _momX_LAB = hcminverse.transform(_momX_HCM);
156        _momY_LAB = hcminverse.transform(_momY_HCM);
157
158        // Find momenta in XCM frame. Note that it is HCM frame that is
159        // boosted, resulting in a sign change later!
160        const bool doXCM = (momX.betaVec().mod2() < 1.);
161        if (doXCM) {
162          const LorentzTransform xcmboost =
163            LorentzTransform::mkFrameTransformFromBeta(momX.betaVec());
164          _momX_XCM = xcmboost.transform(momX);
165          _momY_XCM = xcmboost.transform(momY);
166        }
167
168        for (const Particle& jp : pX) {
169          // Boost from HCM to LAB.
170          FourMomentum lab = hcminverse.transform(jp.momentum());
171          _ePpzX_LAB += lab.E() + dir * lab.pz();
172          _eMpzX_LAB += lab.E() - dir * lab.pz();
173          Particle plab = jp;
174          plab.setMomentum(lab);
175          _pX_LAB.push_back(plab);
176          // Set XCM. Note that since HCM frame is boosted to XCM frame,
177          // we have a sign change
178          if (doXCM) {
179            const LorentzTransform xcmboost =
180              LorentzTransform::mkFrameTransformFromBeta(_momX_HCM.betaVec());
181            FourMomentum xcm = xcmboost.transform(jp.momentum());
182            _ePpzX_XCM += xcm.E() - xcm.pz(); // Sign + => -
183            _eMpzX_XCM += xcm.E() + xcm.pz(); // Sign - => +
184            Particle pxcm = jp;
185            pxcm.setMomentum(xcm);
186            _pX_XCM.push_back(pxcm);
187          }
188        }
189
190        for (const Particle& jp : pY) {
191          // Boost from HCM to LAB
192          FourMomentum lab = hcminverse.transform(jp.momentum());
193          Particle plab = jp;
194          plab.setMomentum(lab);
195          _pY_LAB.push_back(plab);
196          // Boost from HCM to XCM
197          if (doXCM) {
198            const LorentzTransform xcmboost =
199              LorentzTransform::mkFrameTransformFromBeta(_momX_HCM.betaVec());
200            FourMomentum xcm = xcmboost.transform(jp.momentum());
201            Particle pxcm = jp;
202            pxcm.setMomentum(xcm);
203            _pY_XCM.push_back(pxcm);
204          }
205        }
206
207        // Find t: Currently can only handle gap on proton side.
208        // @TODO: Expand to also handle gap on photon side
209        // Boost p from LAB to HCM frame to find t.
210        const FourMomentum proton = hcmboost.transform(diskin.beamHadron().momentum());
211        FourMomentum pPom         = proton - _momY_HCM;
212        _t                        = pPom * pPom;
213
214      }
215
216    private:
217
218      double _M2X, _M2Y, _t;
219      double _gap, _gapUpp, _gapLow;
220      double _ePpzX_LAB, _eMpzX_LAB, _ePpzX_HCM, _eMpzX_HCM, _ePpzX_XCM, _eMpzX_XCM;
221      FourMomentum _momX_HCM, _momY_HCM,_momX_LAB, _momY_LAB, _momX_XCM, _momY_XCM;
222      Particles _pX_HCM, _pY_HCM, _pX_LAB, _pY_LAB, _pX_XCM, _pY_XCM;
223
224    };
225
226    /// Projection to boost system X (photon+Pomeron) particles into its rest frame.
227    ///
228    /// @author Ilkka Helenius
229    class BoostedXSystem : public FinalState {
230    public:
231
232      BoostedXSystem(const FinalState& fs) {
233        setName("BoostedXSystem");
234        declare(fs,"FS");
235        declare(RapidityGap(), "RAPGAP");
236      }
237
238      // Return the boost to XCM frame.
239      const LorentzTransform& boost() const { return _boost; }
240
241      DEFAULT_RIVET_PROJ_CLONE(BoostedXSystem);
242
243    protected:
244
245      // Apply the projection on the supplied event.
246      void project(const Event& e){
247
248        const RapidityGap& rg = apply<RapidityGap>(e, "RAPGAP");
249
250        // Total momentum of the system X.
251        const FourMomentum pX = rg.pX(RapidityGap::HCM);
252
253        // Reset the boost. Is there a separate method for this?
254        _boost = combine(_boost, _boost.inverse());
255
256        // Define boost only when numerically safe, otherwise negligible.
257        if (pX.betaVec().mod2() < 1.)
258          _boost = LorentzTransform::mkFrameTransformFromBeta(pX.betaVec());
259
260        // Boost the particles from system X.
261        _theParticles.clear();
262        _theParticles.reserve(rg.systemX(RapidityGap::HCM).size());
263        for (const Particle& p : rg.systemX(RapidityGap::HCM)) {
264          Particle temp = p;
265          temp.setMomentum(_boost.transform(temp.momentum()));
266          _theParticles.push_back(temp);
267        }
268
269      }
270
271      // Compare projections.
272      CmpState compare(const Projection& p) const {
273        const BoostedXSystem& other = pcast<BoostedXSystem>(p);
274        return mkNamedPCmp(other, "RAPGAP") || mkNamedPCmp(other, "FS");
275      }
276
277    private:
278
279      LorentzTransform _boost;
280
281    };
282
283  }
284
285
286
287  /// @brief H1 diffractive dijets
288  ///
289  /// Diffractive dijets H1 with 920 GeV p and 27.5 GeV e
290  /// Note tagged protons!
291  ///
292  /// @author Christine O. Rasmussen
293  class H1_2007_I746380 : public Analysis {
294  public:
295
296    typedef H1_2007_I746380_PROJECTIONS::RapidityGap RapidityGap;
297    typedef H1_2007_I746380_PROJECTIONS::BoostedXSystem BoostedXSystem;
298
299    /// Constructor
300    RIVET_DEFAULT_ANALYSIS_CTOR(H1_2007_I746380);
301
302    /// @name Analysis methods
303    //@{
304
305    // Book projections and histograms
306    void init() {
307
308      declare(DISKinematics(), "Kinematics");
309      const DISFinalState& disfs = declare(DISFinalState(DISFinalState::BoostFrame::HCM), "DISFS");
310      const BoostedXSystem& disfsXcm = declare( BoostedXSystem(disfs), "BoostedXFS");
311      declare(FastJets(disfsXcm, fastjet::JetAlgorithm::kt_algorithm, fastjet::RecombinationScheme::pt_scheme, 1.0,
312                       JetAlg::Muons::ALL, JetAlg::Invisibles::NONE, nullptr), "DISFSJets");
313      declare(RapidityGap(), "RapidityGap");
314
315      // Book histograms from REF data
316      book(_h_DIS_dsigdzPom, 1, 1, 1);
317      book(_h_DIS_dsigdlogXpom, 2, 1, 1);
318      book(_h_DIS_dsigdW, 3, 1, 1);
319      book(_h_DIS_dsigdQ2, 4, 1, 1);
320      book(_h_DIS_dsigdEtJet1, 5, 1, 1);
321      book(_h_DIS_dsigdAvgEta, 6, 1, 1);
322      book(_h_DIS_dsigdDeltaEta, 7, 1, 1);
323
324      book(_h_PHO_dsigdzPom, 8, 1, 1);
325      book(_h_PHO_dsigdxGam, 9, 1, 1);
326      book(_h_PHO_dsigdlogXpom, 10, 1, 1);
327      book(_h_PHO_dsigdW, 11, 1, 1);
328      book(_h_PHO_dsigdEtJet1, 12, 1, 1);
329      book(_h_PHO_dsigdAvgEta, 13, 1, 1);
330      book(_h_PHO_dsigdDeltaEta, 14, 1, 1);
331      book(_h_PHO_dsigdMjets, 15, 1, 1);
332
333      isDIS  = false;
334      nVeto0 = 0;
335      nVeto1 = 0;
336      nVeto2 = 0;
337      nVeto3 = 0;
338      nVeto4 = 0;
339      nVeto5 = 0;
340      nPHO   = 0;
341      nDIS   = 0;
342    }
343
344    // Do the analysis
345    void analyze(const Event& event) {
346
347      // Event weight
348      isDIS  = false;
349
350      // Projections - special handling of events where no proton found:
351      const RapidityGap&    rg = apply<RapidityGap>(event, "RapidityGap");
352      const DISKinematics& kin = apply<DISKinematics>(event, "Kinematics");
353      const BoostedXSystem& disfsXcm = apply<BoostedXSystem>( event, "BoostedXFS");
354
355      // Determine kinematics: H1 has +z = proton direction
356      int dir   = kin.orientation();
357      double W2 = kin.W2();
358      double W  = sqrt(W2);
359      double y  = kin.y();
360      double Q2 = kin.Q2();
361
362      // Separate into DIS and PHO regimes else veto
363      if (!inRange(W, 165.*GeV, 242.*GeV)) vetoEvent;
364      if (Q2 < 0.01*GeV2) {
365        isDIS = false;
366        ++nPHO;
367      } else if (inRange(Q2, 4.0*GeV2, 80.*GeV2)) {
368        isDIS = true;
369        ++nDIS;
370      } else {
371        vetoEvent;
372      }
373      ++nVeto0;
374
375      // Find diffractive variables as defined in paper.
376      const double M2Y  = rg.M2Y();
377      const double M2X  = rg.M2X();
378      const double abst = abs(rg.t());
379      const double xPom = (isDIS) ? (Q2 + M2X) / (Q2 + W2) :
380        rg.EpPzX(RapidityGap::LAB) / (2. * kin.beamHadron().E());
381
382      // Veto if outside allowed region
383      if (sqrt(M2Y) > 1.6*GeV)    vetoEvent;
384      ++nVeto1;
385      if (abst > 1.0*GeV2) vetoEvent;
386      ++nVeto2;
387      if (xPom > 0.03)     vetoEvent;
388      ++nVeto3;
389
390      // Jet selection. Note jets are found in photon-proton (XCM)
391      // frame, but eta cut is applied in lab frame!
392      Cut jetcuts = Cuts::Et > 4.* GeV;
393      Jets jets   = apply<FastJets>(event, "DISFSJets").jets(jetcuts, cmpMomByEt);
394      // Veto if not dijets and if Et_j1 < 5.0
395      if (jets.size() < 2)       vetoEvent;
396      if (jets[0].Et() < 5.*GeV) vetoEvent;
397      ++nVeto4;
398      // Find Et_jet1 and deltaEta* in XCM frame
399      double EtJet1       = jets[0].Et() * GeV;
400      double etaXCMJet1   = jets[0].eta();
401      double etaXCMJet2   = jets[1].eta();
402      double deltaEtaJets = abs(etaXCMJet1 - etaXCMJet2);
403
404      // Transform from XCM to HCM
405      const LorentzTransform xcmboost = disfsXcm.boost();
406      for (int i = 0; i < 2; ++i) jets[i].transformBy(xcmboost.inverse());
407      // Find mass of jets and EpPz, EmPz of jets
408      FourMomentum momJets = jets[0].momentum() + jets[1].momentum();
409      double M2jets   = momJets.mass2();
410      double EpPzJets = 0.;
411      double EmPzJets = 0.;
412      // DIS variables are found in XCM frame, so boost back again
413      if (isDIS){
414        for (int i = 0; i < 2; ++i) jets[i].transformBy(xcmboost);
415      }
416      // Note sign change wrt. H1 because photon is in +z direction
417      // Jets in HCM so no need to consider orientation.
418      for (int i = 0; i < 2; ++i){
419        EpPzJets += jets[i].E() - jets[i].pz(); // Sign: + => -
420        EmPzJets += jets[i].E() + jets[i].pz(); // Sign: - => +
421      }
422
423      // Transform the jets from HCM to LAB frame where eta cut is
424      // applied for photoproduction.
425      const LorentzTransform hcmboost = kin.boostHCM();
426      for (int i = 0; i < 2; ++i) jets[i].transformBy(hcmboost.inverse());
427      double etaLabJet1 = dir * jets[0].eta();
428      double etaLabJet2 = dir * jets[1].eta();
429      double etaMin     = (isDIS) ? -3. : -1.;
430      double etaMax     = (isDIS) ? 0. : 2.;
431      double eta1       = (isDIS) ? etaXCMJet1 : etaLabJet1;
432      double eta2       = (isDIS) ? etaXCMJet2 : etaLabJet2;
433      if (!inRange(eta1, etaMin, etaMax)) vetoEvent;
434      if (!inRange(eta2, etaMin, etaMax)) vetoEvent;
435      ++nVeto5;
436
437      // Pseudorapidity distributions are examined in lab frame:
438      double avgEtaJets   = 0.5 * (etaLabJet1 + etaLabJet2);
439
440      // Derive xPom and xGam values from the jet kinematics.
441      double zPomJets, xGamJets;
442      if (isDIS) {
443        zPomJets = (Q2 + M2jets) / (Q2 + M2X);
444        xGamJets = EmPzJets / rg.EmPzX(RapidityGap::XCM);
445      } else {
446        // Boost E_p, E_e to HCM frame
447        FourMomentum lep = hcmboost.transform(kin.beamLepton().momentum());
448        FourMomentum had = hcmboost.transform(kin.beamHadron().momentum());
449        zPomJets = EpPzJets / (2. * xPom * had.E());
450        xGamJets = EmPzJets / (2. * y * lep.E());
451      }
452
453      // Now fill histograms
454      if (isDIS){
455        _h_DIS_dsigdzPom     ->fill(zPomJets);
456        _h_DIS_dsigdlogXpom  ->fill(log10(xPom));
457        _h_DIS_dsigdW        ->fill(W);
458        _h_DIS_dsigdQ2       ->fill(Q2);
459        _h_DIS_dsigdEtJet1   ->fill(EtJet1);
460        _h_DIS_dsigdAvgEta   ->fill(avgEtaJets);
461        _h_DIS_dsigdDeltaEta ->fill(deltaEtaJets);
462      } else {
463        _h_PHO_dsigdzPom     ->fill(zPomJets);
464        _h_PHO_dsigdxGam     ->fill(xGamJets);
465        _h_PHO_dsigdlogXpom  ->fill(log10(xPom));
466        _h_PHO_dsigdW        ->fill(W);
467        _h_PHO_dsigdEtJet1   ->fill(EtJet1);
468        _h_PHO_dsigdAvgEta   ->fill(avgEtaJets);
469        _h_PHO_dsigdDeltaEta ->fill(deltaEtaJets);
470        _h_PHO_dsigdMjets    ->fill(sqrt(M2jets));
471      }
472
473    }
474
475    // Finalize
476    void finalize() {
477      // Normalise to cross section
478      const double norm = crossSection()/picobarn/sumOfWeights();
479
480      scale( _h_DIS_dsigdzPom    , norm);
481      scale( _h_DIS_dsigdlogXpom , norm);
482      scale( _h_DIS_dsigdW       , norm);
483      scale( _h_DIS_dsigdQ2      , norm);
484      scale( _h_DIS_dsigdEtJet1  , norm);
485      scale( _h_DIS_dsigdAvgEta  , norm);
486      scale( _h_DIS_dsigdDeltaEta, norm);
487
488      scale( _h_PHO_dsigdzPom    , norm);
489      scale( _h_PHO_dsigdxGam    , norm);
490      scale( _h_PHO_dsigdlogXpom , norm);
491      scale( _h_PHO_dsigdW       , norm);
492      scale( _h_PHO_dsigdEtJet1  , norm);
493      scale( _h_PHO_dsigdAvgEta  , norm);
494      scale( _h_PHO_dsigdDeltaEta, norm);
495      scale( _h_PHO_dsigdMjets   , norm);
496
497      const double dPHO = nPHO;
498      MSG_INFO("H1_2007_I746380");
499      MSG_INFO("Cross section = " << crossSection()/picobarn << " pb");
500      MSG_INFO("Number of events = " << numEvents() << ", sumW = " << sumOfWeights());
501      MSG_INFO("Number of PHO = " << nPHO << ", number of DIS = " << nDIS);
502      MSG_INFO("Events passing electron veto   = " << nVeto0 << " (" << nVeto0/dPHO * 100. << "%)" );
503      MSG_INFO("Events passing MY              = " << nVeto1 << " (" << nVeto1/dPHO * 100. << "%)" );
504      MSG_INFO("Events passing t veto          = " << nVeto2 << " (" << nVeto2/dPHO * 100. << "%)" );
505      MSG_INFO("Events passing xPom            = " << nVeto3 << " (" << nVeto3/dPHO * 100. << "%)" );
506      MSG_INFO("Events passing jet Et veto     = " << nVeto4 << " (" << nVeto4/dPHO * 100. << "%)" );
507      MSG_INFO("Events passing jet eta veto    = " << nVeto5 << " (" << nVeto5/dPHO * 100. << "%)" );
508
509    }
510
511    //@}
512
513
514  private:
515
516    /// @name Histograms
517    //@{
518    // Book histograms from REF data
519    Histo1DPtr _h_DIS_dsigdzPom    ;
520    Histo1DPtr _h_DIS_dsigdlogXpom ;
521    Histo1DPtr _h_DIS_dsigdW       ;
522    Histo1DPtr _h_DIS_dsigdQ2      ;
523    Histo1DPtr _h_DIS_dsigdEtJet1  ;
524    Histo1DPtr _h_DIS_dsigdAvgEta  ;
525    Histo1DPtr _h_DIS_dsigdDeltaEta;
526
527    Histo1DPtr _h_PHO_dsigdzPom    ;
528    Histo1DPtr _h_PHO_dsigdxGam    ;
529    Histo1DPtr _h_PHO_dsigdlogXpom ;
530    Histo1DPtr _h_PHO_dsigdW       ;
531    Histo1DPtr _h_PHO_dsigdEtJet1  ;
532    Histo1DPtr _h_PHO_dsigdAvgEta  ;
533    Histo1DPtr _h_PHO_dsigdDeltaEta;
534    Histo1DPtr _h_PHO_dsigdMjets   ;
535    //@}
536
537    bool isDIS;
538    int  nVeto0, nVeto1, nVeto2, nVeto3, nVeto4, nVeto5;
539    int nPHO, nDIS;
540  };
541
542  RIVET_DECLARE_PLUGIN(H1_2007_I746380);
543
544}