1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/UnstableParticles.hh"
  4
  5namespace Rivet {
  6
  7
  8  /// @brief chi_c production at 7 TeV
  9  class ATLAS_2014_I1292798 : public Analysis {
 10  public:
 11
 12    /// Constructor
 13    RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2014_I1292798);
 14
 15
 16    /// @name Analysis methods
 17    /// @{
 18
 19    /// Book histograms and initialise projections before the run
 20    void init() {
 21      declare(UnstableParticles(Cuts::pid==443 || Cuts::pid==20443 || Cuts::pid==445), "UFS");
 22      for (unsigned int ichi=0; ichi<2; ++ichi) {
 23        for (unsigned int ipT=0; ipT<2; ++ipT) {
 24          for (unsigned int iprompt=0; iprompt<2; ++iprompt) {
 25            book(_h_pT[ichi][iprompt][ipT], 1+ichi+2*iprompt+4*ipT, 1, 1);
 26          }
 27          book(_h_pT[ichi][2][ipT],
 28               "TMP/chi_"+toString(ichi)+"_"+toString(ipT),
 29               refData(1+ichi+4*ipT,1,1));
 30        }
 31      }
 32      for (unsigned int iprompt=0; iprompt<2; ++iprompt) {
 33        for (unsigned int ifeed=0; ifeed<2; ++ifeed) {
 34          book(_h_Jpsi[iprompt][ifeed],
 35               "TMP/jpsi_"+toString(iprompt)+"_"+toString(ifeed),
 36               refData(9+iprompt,1,1));
 37        }
 38      }
 39    }
 40
 41
 42    /// Perform the per-event analysis
 43    void analyze(const Event& event) {
 44
 45      // Final state of unstable particles to get particle spectra
 46      const UnstableParticles& ufs = apply<UnstableParticles>(event, "UFS");
 47
 48      for (const Particle& p : ufs.particles()) {
 49        // prompt/non-prompt
 50        bool nonPrompt = p.fromBottom();
 51        // J/psi as a reference
 52        if (p.pid()==443) {
 53          if (p.absrap()<0.75) _h_Jpsi[nonPrompt][0]->fill(p.perp());
 54        }
 55        else if (p.children().size()==2) {
 56          Particle Jpsi;
 57          if(p.children()[0].pid()==22 && p.children()[1].pid()==443) {
 58            Jpsi=p.children()[1];
 59          }
 60          else if(p.children()[1].pid()==22 && p.children()[0].pid()==443) {
 61            Jpsi=p.children()[0];
 62          }
 63          else  continue;
 64          if (Jpsi.absrap()>0.75) continue;
 65          unsigned int ichi = p.pid()==20443 ? 0 : 1;
 66          _h_Jpsi[nonPrompt][1]->fill(Jpsi.perp());
 67          _h_pT[ichi][nonPrompt][0]->fill(Jpsi.perp());
 68          _h_pT[ichi][nonPrompt][1]->fill(p   .perp());
 69          _h_pT[ichi][2        ][0]->fill(Jpsi.perp());
 70          _h_pT[ichi][2        ][1]->fill(p   .perp());
 71        }
 72      }
 73    }
 74
 75
 76    /// Normalise histograms etc., after the run
 77    void finalize() {
 78      // J/psi -> mu+mu- branching ratio
 79      const double br = 0.05961;
 80      const double fact = br*crossSection()/picobarn/sumOfWeights();
 81      for (unsigned int iprompt=0; iprompt<3; ++iprompt) {
 82        for (unsigned int ichi=0; ichi<2; ++ichi) {
 83          for (unsigned int ipT=0; ipT<2; ++ipT) {
 84            scale(_h_pT[ichi][iprompt][ipT],fact);
 85          }
 86        }
 87      }
 88      for (unsigned int iprompt=0; iprompt<2; ++iprompt) {
 89        for (unsigned int ifeed=0; ifeed<2; ++ifeed) {
 90          scale(_h_Jpsi[iprompt][ifeed], fact);
 91        }
 92      }
 93      // Jpsi feed down prompt
 94      Estimate1DPtr tmp;
 95      book(tmp,9,1,1);
 96      efficiency(_h_Jpsi[0][1], _h_Jpsi[0][0], tmp);
 97      // ratio prompt chi2/chi1
 98      book(tmp,10,1,1);
 99      divide(_h_pT[1][0][0], _h_pT[0][0][0], tmp);
100      // ratio non-prompt chi2/chi1
101      book(tmp,11,1,1);
102      divide(_h_pT[1][1][0], _h_pT[0][1][0], tmp);
103      // ratio non-prmopt chi1
104      book(tmp,12,1,1);
105      divide(_h_pT[0][1][1], _h_pT[0][2][1], tmp);
106      // ratio non-prmopt chi2
107      book(tmp,13,1,1);
108      divide(_h_pT[1][1][1], _h_pT[1][2][1], tmp);
109    }
110
111    /// @}
112
113
114    /// @name Histograms
115    /// @{
116    Histo1DPtr _h_pT[2][3][2],_h_Jpsi[2][2];
117    /// @}
118
119
120  };
121
122
123  RIVET_DECLARE_PLUGIN(ATLAS_2014_I1292798);
124
125}