Page MenuHomeHEPForge

No OneTemporary

diff --git a/analyses/pluginATLAS/ATLAS_2011_I945498.cc b/analyses/pluginATLAS/ATLAS_2011_I945498.cc
--- a/analyses/pluginATLAS/ATLAS_2011_I945498.cc
+++ b/analyses/pluginATLAS/ATLAS_2011_I945498.cc
@@ -1,305 +1,303 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/ZFinder.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/VetoedFinalState.hh"
#include "Rivet/Projections/IdentifiedFinalState.hh"
#include "Rivet/Projections/LeadingParticlesFinalState.hh"
namespace Rivet {
/// ATLAS Z+jets in pp at 7 TeV
class ATLAS_2011_I945498 : public Analysis {
public:
/// Constructor
ATLAS_2011_I945498()
: Analysis("ATLAS_2011_I945498")
{ }
/// Book histograms and initialise projections before the run
void init() {
// Variable initialisation
_isZeeSample = false;
_isZmmSample = false;
for (size_t chn = 0; chn < 3; ++chn) {
- weights_nj0[chn] = 0;
- weights_nj1[chn] = 0;
- weights_nj2[chn] = 0;
- weights_nj3[chn] = 0;
- weights_nj4[chn] = 0;
+ book(weights_nj0[chn], (ostringstream("weights_nj0_") << chn).str());
+ book(weights_nj1[chn], (ostringstream("weights_nj1_") << chn).str());
+ book(weights_nj2[chn], (ostringstream("weights_nj2_") << chn).str());
+ book(weights_nj3[chn], (ostringstream("weights_nj3_") << chn).str());
+ book(weights_nj4[chn], (ostringstream("weights_nj4_") << chn).str());
}
// Set up projections
FinalState fs;
ZFinder zfinder_mu(fs, Cuts::abseta < 2.4 && Cuts::pT > 20*GeV, PID::MUON, 66*GeV, 116*GeV, 0.1, ZFinder::CLUSTERNODECAY);
declare(zfinder_mu, "ZFinder_mu");
Cut cuts = (Cuts::abseta < 1.37 || Cuts::absetaIn(1.52, 2.47)) && Cuts::pT > 20*GeV;
ZFinder zfinder_el(fs, cuts, PID::ELECTRON, 66*GeV, 116*GeV, 0.1, ZFinder::CLUSTERNODECAY);
declare(zfinder_el, "ZFinder_el");
Cut cuts25_20 = Cuts::abseta < 2.5 && Cuts::pT > 20*GeV;
// For combined cross-sections (combined phase space + dressed level)
ZFinder zfinder_comb_mu(fs, cuts25_20, PID::MUON, 66.0*GeV, 116.0*GeV, 0.1, ZFinder::CLUSTERNODECAY);
declare(zfinder_comb_mu, "ZFinder_comb_mu");
ZFinder zfinder_comb_el(fs, cuts25_20, PID::ELECTRON, 66.0*GeV, 116.0*GeV, 0.1, ZFinder::CLUSTERNODECAY);
declare(zfinder_comb_el, "ZFinder_comb_el");
// Define veto FS in order to prevent Z-decay products entering the jet algorithm
VetoedFinalState remfs;
remfs.addVetoOnThisFinalState(zfinder_el);
remfs.addVetoOnThisFinalState(zfinder_mu);
VetoedFinalState remfs_comb;
remfs_comb.addVetoOnThisFinalState(zfinder_comb_el);
remfs_comb.addVetoOnThisFinalState(zfinder_comb_mu);
FastJets jets(remfs, FastJets::ANTIKT, 0.4);
jets.useInvisibles();
declare(jets, "jets");
FastJets jets_comb(remfs_comb, FastJets::ANTIKT, 0.4);
jets_comb.useInvisibles();
declare(jets_comb, "jets_comb");
// 0=el, 1=mu, 2=comb
for (size_t chn = 0; chn < 3; ++chn) {
book(_h_njet_incl[chn] ,1, 1, chn+1);
book(_h_njet_ratio[chn] ,2, 1, chn+1);
book(_h_ptjet[chn] ,3, 1, chn+1);
book(_h_ptlead[chn] ,4, 1, chn+1);
book(_h_ptseclead[chn] ,5, 1, chn+1);
book(_h_yjet[chn] ,6, 1, chn+1);
book(_h_ylead[chn] ,7, 1, chn+1);
book(_h_yseclead[chn] ,8, 1, chn+1);
book(_h_mass[chn] ,9, 1, chn+1);
book(_h_deltay[chn] ,10, 1, chn+1);
book(_h_deltaphi[chn] ,11, 1, chn+1);
book(_h_deltaR[chn] ,12, 1, chn+1);
}
}
// Jet selection criteria universal for electron and muon channel
/// @todo Replace with a Cut passed to jetsByPt
Jets selectJets(const ZFinder* zf, const FastJets* allJets) {
const FourMomentum l1 = zf->constituents()[0].momentum();
const FourMomentum l2 = zf->constituents()[1].momentum();
Jets jets;
foreach (const Jet& jet, allJets->jetsByPt(30*GeV)) {
const FourMomentum jmom = jet.momentum();
if (jmom.absrap() < 4.4 &&
deltaR(l1, jmom) > 0.5 && deltaR(l2, jmom) > 0.5) {
jets.push_back(jet);
}
}
return jets;
}
/// Perform the per-event analysis
void analyze(const Event& event) {
- const double weight = 1.0;
-
vector<const ZFinder*> zfs;
zfs.push_back(& (apply<ZFinder>(event, "ZFinder_el")));
zfs.push_back(& (apply<ZFinder>(event, "ZFinder_mu")));
zfs.push_back(& (apply<ZFinder>(event, "ZFinder_comb_el")));
zfs.push_back(& (apply<ZFinder>(event, "ZFinder_comb_mu")));
vector<const FastJets*> fjs;
fjs.push_back(& (apply<FastJets>(event, "jets")));
fjs.push_back(& (apply<FastJets>(event, "jets_comb")));
// Determine what kind of MC sample this is
const bool isZee = (zfs[0]->bosons().size() == 1) || (zfs[2]->bosons().size() == 1);
const bool isZmm = (zfs[1]->bosons().size() == 1) || (zfs[3]->bosons().size() == 1);
if (isZee) _isZeeSample = true;
if (isZmm) _isZmmSample = true;
// Require exactly one electronic or muonic Z-decay in the event
bool isZeemm = ( (zfs[0]->bosons().size() == 1 && zfs[1]->bosons().size() != 1) ||
(zfs[1]->bosons().size() == 1 && zfs[0]->bosons().size() != 1) );
bool isZcomb = ( (zfs[2]->bosons().size() == 1 && zfs[3]->bosons().size() != 1) ||
(zfs[3]->bosons().size() == 1 && zfs[2]->bosons().size() != 1) );
if (!isZeemm && !isZcomb) vetoEvent;
vector<int> zfIDs;
vector<int> fjIDs;
if (isZeemm) {
int chn = zfs[0]->bosons().size() == 1 ? 0 : 1;
zfIDs.push_back(chn);
fjIDs.push_back(0);
}
if (isZcomb) {
int chn = zfs[2]->bosons().size() == 1 ? 2 : 3;
zfIDs.push_back(chn);
fjIDs.push_back(1);
}
for (size_t izf = 0; izf < zfIDs.size(); ++izf) {
int zfID = zfIDs[izf];
int fjID = fjIDs[izf];
int chn = zfID;
if (zfID == 2 || zfID == 3) chn = 2;
Jets jets = selectJets(zfs[zfID], fjs[fjID]);
switch (jets.size()) {
case 0:
- weights_nj0[chn] += weight;
+ weights_nj0[chn]->fill();
break;
case 1:
- weights_nj0[chn] += weight;
- weights_nj1[chn] += weight;
+ weights_nj0[chn]->fill();
+ weights_nj1[chn]->fill();
break;
case 2:
- weights_nj0[chn] += weight;
- weights_nj1[chn] += weight;
- weights_nj2[chn] += weight;
+ weights_nj0[chn]->fill();
+ weights_nj1[chn]->fill();
+ weights_nj2[chn]->fill();
break;
case 3:
- weights_nj0[chn] += weight;
- weights_nj1[chn] += weight;
- weights_nj2[chn] += weight;
- weights_nj3[chn] += weight;
+ weights_nj0[chn]->fill();
+ weights_nj1[chn]->fill();
+ weights_nj2[chn]->fill();
+ weights_nj3[chn]->fill();
break;
default: // >= 4
- weights_nj0[chn] += weight;
- weights_nj1[chn] += weight;
- weights_nj2[chn] += weight;
- weights_nj3[chn] += weight;
- weights_nj4[chn] += weight;
+ weights_nj0[chn]->fill();
+ weights_nj1[chn]->fill();
+ weights_nj2[chn]->fill();
+ weights_nj3[chn]->fill();
+ weights_nj4[chn]->fill();
}
// Require at least one jet
if (jets.empty()) continue;
// Fill jet multiplicities
for (size_t ijet = 1; ijet <= jets.size(); ++ijet) {
- _h_njet_incl[chn]->fill(ijet, weight);
+ _h_njet_incl[chn]->fill(ijet);
}
// Loop over selected jets, fill inclusive jet distributions
for (size_t ijet = 0; ijet < jets.size(); ++ijet) {
- _h_ptjet[chn]->fill(jets[ijet].pT()/GeV, weight);
- _h_yjet [chn]->fill(fabs(jets[ijet].rapidity()), weight);
+ _h_ptjet[chn]->fill(jets[ijet].pT()/GeV);
+ _h_yjet [chn]->fill(fabs(jets[ijet].rapidity()));
}
// Leading jet histos
const double ptlead = jets[0].pT()/GeV;
const double yabslead = fabs(jets[0].rapidity());
- _h_ptlead[chn]->fill(ptlead, weight);
- _h_ylead [chn]->fill(yabslead, weight);
+ _h_ptlead[chn]->fill(ptlead);
+ _h_ylead [chn]->fill(yabslead);
if (jets.size() >= 2) {
// Second jet histos
const double pt2ndlead = jets[1].pT()/GeV;
const double yabs2ndlead = fabs(jets[1].rapidity());
- _h_ptseclead[chn] ->fill(pt2ndlead, weight);
- _h_yseclead [chn] ->fill(yabs2ndlead, weight);
+ _h_ptseclead[chn] ->fill(pt2ndlead);
+ _h_yseclead [chn] ->fill(yabs2ndlead);
// Dijet histos
const double deltaphi = fabs(deltaPhi(jets[1], jets[0]));
const double deltarap = fabs(jets[0].rapidity() - jets[1].rapidity()) ;
const double deltar = fabs(deltaR(jets[0], jets[1], RAPIDITY));
const double mass = (jets[0].momentum() + jets[1].momentum()).mass();
- _h_mass [chn] ->fill(mass/GeV, weight);
- _h_deltay [chn] ->fill(deltarap, weight);
- _h_deltaphi[chn] ->fill(deltaphi, weight);
- _h_deltaR [chn] ->fill(deltar, weight);
+ _h_mass [chn] ->fill(mass/GeV);
+ _h_deltay [chn] ->fill(deltarap);
+ _h_deltaphi[chn] ->fill(deltaphi);
+ _h_deltaR [chn] ->fill(deltar);
}
}
}
/// @name Ratio calculator util functions
//@{
/// Calculate the ratio, being careful about div-by-zero
double ratio(double a, double b) {
return (b != 0) ? a/b : 0;
}
/// Calculate the ratio error, being careful about div-by-zero
double ratio_err(double a, double b) {
return (b != 0) ? sqrt(a/sqr(b) + sqr(a)/(b*b*b)) : 0;
}
//@}
void finalize() {
// Fill ratio histograms
for (size_t chn = 0; chn < 3; ++chn) {
- _h_njet_ratio[chn]->addPoint(1, ratio(weights_nj1[chn], weights_nj0[chn]), 0.5, ratio_err(weights_nj1[chn], weights_nj0[chn]));
- _h_njet_ratio[chn]->addPoint(2, ratio(weights_nj2[chn], weights_nj1[chn]), 0.5, ratio_err(weights_nj2[chn], weights_nj1[chn]));
- _h_njet_ratio[chn]->addPoint(3, ratio(weights_nj3[chn], weights_nj2[chn]), 0.5, ratio_err(weights_nj3[chn], weights_nj2[chn]));
- _h_njet_ratio[chn]->addPoint(4, ratio(weights_nj4[chn], weights_nj3[chn]), 0.5, ratio_err(weights_nj4[chn], weights_nj3[chn]));
+ _h_njet_ratio[chn]->addPoint(1, ratio(weights_nj1[chn]->val(), weights_nj0[chn]->val()), 0.5, ratio_err(weights_nj1[chn]->val(), weights_nj0[chn]->val()));
+ _h_njet_ratio[chn]->addPoint(2, ratio(weights_nj2[chn]->val(), weights_nj1[chn]->val()), 0.5, ratio_err(weights_nj2[chn]->val(), weights_nj1[chn]->val()));
+ _h_njet_ratio[chn]->addPoint(3, ratio(weights_nj3[chn]->val(), weights_nj2[chn]->val()), 0.5, ratio_err(weights_nj3[chn]->val(), weights_nj2[chn]->val()));
+ _h_njet_ratio[chn]->addPoint(4, ratio(weights_nj4[chn]->val(), weights_nj3[chn]->val()), 0.5, ratio_err(weights_nj4[chn]->val(), weights_nj3[chn]->val()));
}
// Scale other histos
for (size_t chn = 0; chn < 3; ++chn) {
// For ee and mumu channels: normalize to Njet inclusive cross-section
double xs = (chn == 2) ? crossSectionPerEvent()/picobarn : 1 / weights_nj0[chn];
// For inclusive MC sample(ee/mmu channels together) we want the single-lepton-flavor xsec
if (_isZeeSample && _isZmmSample) xs /= 2;
// Special case histogram: always not normalized
scale(_h_njet_incl[chn], (chn < 2) ? crossSectionPerEvent()/picobarn : xs);
scale(_h_ptjet[chn] , xs);
scale(_h_ptlead[chn] , xs);
scale(_h_ptseclead[chn], xs);
scale(_h_yjet[chn] , xs);
scale(_h_ylead[chn] , xs);
scale(_h_yseclead[chn] , xs);
scale(_h_deltaphi[chn] , xs);
scale(_h_deltay[chn] , xs);
scale(_h_deltaR[chn] , xs);
scale(_h_mass[chn] , xs);
}
}
//@}
private:
bool _isZeeSample;
bool _isZmmSample;
- double weights_nj0[3];
- double weights_nj1[3];
- double weights_nj2[3];
- double weights_nj3[3];
- double weights_nj4[3];
+ CounterPtr weights_nj0[3];
+ CounterPtr weights_nj1[3];
+ CounterPtr weights_nj2[3];
+ CounterPtr weights_nj3[3];
+ CounterPtr weights_nj4[3];
Scatter2DPtr _h_njet_ratio[3];
Histo1DPtr _h_njet_incl[3];
Histo1DPtr _h_ptjet[3];
Histo1DPtr _h_ptlead[3];
Histo1DPtr _h_ptseclead[3];
Histo1DPtr _h_yjet[3];
Histo1DPtr _h_ylead[3];
Histo1DPtr _h_yseclead[3];
Histo1DPtr _h_deltaphi[3];
Histo1DPtr _h_deltay[3];
Histo1DPtr _h_deltaR[3];
Histo1DPtr _h_mass[3];
};
DECLARE_RIVET_PLUGIN(ATLAS_2011_I945498);
}
diff --git a/analyses/pluginATLAS/ATLAS_2011_S8924791.cc b/analyses/pluginATLAS/ATLAS_2011_S8924791.cc
--- a/analyses/pluginATLAS/ATLAS_2011_S8924791.cc
+++ b/analyses/pluginATLAS/ATLAS_2011_S8924791.cc
@@ -1,127 +1,125 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Projections/VetoedFinalState.hh"
#include "Rivet/Projections/VisibleFinalState.hh"
#include "Rivet/Projections/JetShape.hh"
namespace Rivet {
/// @brief ATLAS jet shape analysis
/// @author Andy Buckley, Judith Katzy, Francesc Vives
class ATLAS_2011_S8924791 : public Analysis {
public:
/// Constructor
ATLAS_2011_S8924791()
: Analysis("ATLAS_2011_S8924791")
{ }
/// @name Analysis methods
//@{
void init() {
// Set up projections
const FinalState fs(-5.0, 5.0);
declare(fs, "FS");
FastJets fj(fs, FastJets::ANTIKT, 0.6);
fj.useInvisibles();
declare(fj, "Jets");
// Specify pT bins
_ptedges = {{ 30.0, 40.0, 60.0, 80.0, 110.0, 160.0, 210.0, 260.0, 310.0, 400.0, 500.0, 600.0 }};
_yedges = {{ 0.0, 0.3, 0.8, 1.2, 2.1, 2.8 }};
// Register a jet shape projection and histogram for each pT bin
for (size_t i = 0; i < 11; ++i) {
for (size_t j = 0; j < 6; ++j) {
if (i == 8 && j == 4) continue;
if (i == 9 && j == 4) continue;
if (i == 10 && j != 5) continue;
// Set up projections for each (pT,y) bin
_jsnames_pT[i][j] = "JetShape" + to_str(i) + to_str(j);
const double ylow = (j < 5) ? _yedges[j] : _yedges.front();
const double yhigh = (j < 5) ? _yedges[j+1] : _yedges.back();
const JetShape jsp(fj, 0.0, 0.7, 7, _ptedges[i], _ptedges[i+1], ylow, yhigh, RAPIDITY);
declare(jsp, _jsnames_pT[i][j]);
// Book profile histograms for each (pT,y) bin
book(_profhistRho_pT[i][j] ,i+1, j+1, 1);
book(_profhistPsi_pT[i][j] ,i+1, j+1, 2);
}
}
}
/// Do the analysis
void analyze(const Event& evt) {
// Get jets and require at least one to pass pT and y cuts
const Jets jets = apply<FastJets>(evt, "Jets")
.jetsByPt(Cuts::ptIn(_ptedges.front()*GeV, _ptedges.back()*GeV) && Cuts::absrap < 2.8);
MSG_DEBUG("Jet multiplicity before cuts = " << jets.size());
if (jets.size() == 0) {
MSG_DEBUG("No jets found in required pT and rapidity range");
vetoEvent;
}
// Calculate and histogram jet shapes
- const double weight = 1.0;
-
for (size_t ipt = 0; ipt < 11; ++ipt) {
for (size_t jy = 0; jy < 6; ++jy) {
if (ipt == 8 && jy == 4) continue;
if (ipt == 9 && jy == 4) continue;
if (ipt == 10 && jy != 5) continue;
const JetShape& jsipt = apply<JetShape>(evt, _jsnames_pT[ipt][jy]);
for (size_t ijet = 0; ijet < jsipt.numJets(); ++ijet) {
for (size_t rbin = 0; rbin < jsipt.numBins(); ++rbin) {
const double r_rho = jsipt.rBinMid(rbin);
- _profhistRho_pT[ipt][jy]->fill(r_rho, (1./0.1)*jsipt.diffJetShape(ijet, rbin), weight);
+ _profhistRho_pT[ipt][jy]->fill(r_rho, (1./0.1)*jsipt.diffJetShape(ijet, rbin));
const double r_Psi = jsipt.rBinMid(rbin);
- _profhistPsi_pT[ipt][jy]->fill(r_Psi, jsipt.intJetShape(ijet, rbin), weight);
+ _profhistPsi_pT[ipt][jy]->fill(r_Psi, jsipt.intJetShape(ijet, rbin));
}
}
}
}
}
// Finalize
void finalize() {
}
//@}
private:
/// @name Analysis data
//@{
/// Jet \f$ p_\perp\f$ bins.
vector<double> _ptedges; // This can't be a raw array if we want to initialise it non-painfully
vector<double> _yedges;
/// JetShape projection name for each \f$p_\perp\f$ bin.
string _jsnames_pT[11][6];
//@}
/// @name Histograms
//@{
Profile1DPtr _profhistRho_pT[11][6];
Profile1DPtr _profhistPsi_pT[11][6];
//@}
};
// The hook for the plugin system
DECLARE_RIVET_PLUGIN(ATLAS_2011_S8924791);
}
diff --git a/analyses/pluginATLAS/ATLAS_2011_S8994773.cc b/analyses/pluginATLAS/ATLAS_2011_S8994773.cc
--- a/analyses/pluginATLAS/ATLAS_2011_S8994773.cc
+++ b/analyses/pluginATLAS/ATLAS_2011_S8994773.cc
@@ -1,138 +1,136 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"
namespace Rivet {
/// @author Jinlong Zhang
class ATLAS_2011_S8994773 : public Analysis {
public:
ATLAS_2011_S8994773()
: Analysis("ATLAS_2011_S8994773") { }
void init() {
const FinalState fs500(-2.5, 2.5, 500*MeV);
declare(fs500, "FS500");
const FinalState fslead(-2.5, 2.5, 1.0*GeV);
declare(fslead, "FSlead");
// Get an index for the beam energy
isqrts = -1;
if (fuzzyEquals(sqrtS(), 900*GeV)) isqrts = 0;
else if (fuzzyEquals(sqrtS(), 7*TeV)) isqrts = 1;
assert(isqrts >= 0);
// N profiles, 500 MeV pT cut
book(_hist_N_transverse_500 ,1+isqrts, 1, 1);
// pTsum profiles, 500 MeV pT cut
book(_hist_ptsum_transverse_500 ,3+isqrts, 1, 1);
// N vs. Delta(phi) profiles, 500 MeV pT cut
book(_hist_N_vs_dPhi_1_500 ,13+isqrts, 1, 1);
book(_hist_N_vs_dPhi_2_500 ,13+isqrts, 1, 2);
book(_hist_N_vs_dPhi_3_500 ,13+isqrts, 1, 3);
}
void analyze(const Event& event) {
- const double weight = 1.0;
-
// Require at least one cluster in the event with pT >= 1 GeV
const FinalState& fslead = apply<FinalState>(event, "FSlead");
if (fslead.size() < 1) {
vetoEvent;
}
// These are the particles with pT > 500 MeV
const FinalState& chargedNeutral500 = apply<FinalState>(event, "FS500");
// Identify leading object and its phi and pT
Particles particles500 = chargedNeutral500.particlesByPt();
Particle p_lead = particles500[0];
const double philead = p_lead.phi();
const double etalead = p_lead.eta();
const double pTlead = p_lead.pT();
MSG_DEBUG("Leading object: pT = " << pTlead << ", eta = " << etalead << ", phi = " << philead);
// Iterate over all > 500 MeV particles and count particles and scalar pTsum in the three regions
vector<double> num500(3, 0), ptSum500(3, 0.0);
// Temporary histos that bin N in dPhi.
// NB. Only one of each needed since binnings are the same for the energies and pT cuts
Histo1D hist_num_dphi_500(refData(13+isqrts,1,1));
foreach (const Particle& p, particles500) {
const double pT = p.pT();
const double dPhi = deltaPhi(philead, p.phi());
const int ir = region_index(dPhi);
num500[ir] += 1;
ptSum500[ir] += pT;
// Fill temp histos to bin N in dPhi
if (p.genParticle() != p_lead.genParticle()) { // We don't want to fill all those zeros from the leading track...
hist_num_dphi_500.fill(dPhi, 1);
}
}
// Now fill underlying event histograms
// The densities are calculated by dividing the UE properties by dEta*dPhi
// -- each region has a dPhi of 2*PI/3 and dEta is two times 2.5
const double dEtadPhi = (2*2.5 * 2*PI/3.0);
- _hist_N_transverse_500->fill(pTlead/GeV, num500[1]/dEtadPhi, weight);
- _hist_ptsum_transverse_500->fill(pTlead/GeV, ptSum500[1]/GeV/dEtadPhi, weight);
+ _hist_N_transverse_500->fill(pTlead/GeV, num500[1]/dEtadPhi);
+ _hist_ptsum_transverse_500->fill(pTlead/GeV, ptSum500[1]/GeV/dEtadPhi);
// Update the "proper" dphi profile histograms
// Note that we fill dN/dEtadPhi: dEta = 2*2.5, dPhi = 2*PI/nBins
// The values tabulated in the note are for an (undefined) signed Delta(phi) rather than
// |Delta(phi)| and so differ by a factor of 2: we have to actually norm for angular range = 2pi
const size_t nbins = refData(13+isqrts,1,1).numPoints();
for (size_t i = 0; i < nbins; ++i) {
double mean = hist_num_dphi_500.bin(i).xMid();
double value = 0.;
if (hist_num_dphi_500.bin(i).numEntries() > 0) {
mean = hist_num_dphi_500.bin(i).xMean();
value = hist_num_dphi_500.bin(i).area()/hist_num_dphi_500.bin(i).xWidth()/10.0;
}
- if (pTlead/GeV >= 1.0) _hist_N_vs_dPhi_1_500->fill(mean, value, weight);
- if (pTlead/GeV >= 2.0) _hist_N_vs_dPhi_2_500->fill(mean, value, weight);
- if (pTlead/GeV >= 3.0) _hist_N_vs_dPhi_3_500->fill(mean, value, weight);
+ if (pTlead/GeV >= 1.0) _hist_N_vs_dPhi_1_500->fill(mean, value);
+ if (pTlead/GeV >= 2.0) _hist_N_vs_dPhi_2_500->fill(mean, value);
+ if (pTlead/GeV >= 3.0) _hist_N_vs_dPhi_3_500->fill(mean, value);
}
}
void finalize() {
}
private:
// Little helper function to identify Delta(phi) regions
inline int region_index(double dphi) {
assert(inRange(dphi, 0.0, PI, CLOSED, CLOSED));
if (dphi < PI/3.0) return 0;
if (dphi < 2*PI/3.0) return 1;
return 2;
}
private:
int isqrts;
Profile1DPtr _hist_N_transverse_500;
Profile1DPtr _hist_ptsum_transverse_500;
Profile1DPtr _hist_N_vs_dPhi_1_500;
Profile1DPtr _hist_N_vs_dPhi_2_500;
Profile1DPtr _hist_N_vs_dPhi_3_500;
};
// The hook for the plugin system
DECLARE_RIVET_PLUGIN(ATLAS_2011_S8994773);
}
diff --git a/analyses/pluginATLAS/ATLAS_2011_S9035664.cc b/analyses/pluginATLAS/ATLAS_2011_S9035664.cc
--- a/analyses/pluginATLAS/ATLAS_2011_S9035664.cc
+++ b/analyses/pluginATLAS/ATLAS_2011_S9035664.cc
@@ -1,138 +1,133 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/Beam.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/ChargedFinalState.hh"
#include "Rivet/Projections/UnstableFinalState.hh"
namespace Rivet {
/// @brief J/psi production at ATLAS
class ATLAS_2011_S9035664: public Analysis {
public:
/// Constructor
ATLAS_2011_S9035664()
: Analysis("ATLAS_2011_S9035664")
{}
/// @name Analysis methods
//@{
void init() {
declare(UnstableFinalState(), "UFS");
book(_nonPrRapHigh , 14, 1, 1);
book(_nonPrRapMedHigh , 13, 1, 1);
book(_nonPrRapMedLow , 12, 1, 1);
book(_nonPrRapLow , 11, 1, 1);
book(_PrRapHigh , 18, 1, 1);
book(_PrRapMedHigh , 17, 1, 1);
book(_PrRapMedLow , 16, 1, 1);
book(_PrRapLow , 15, 1, 1);
book(_IncRapHigh , 20, 1, 1);
book(_IncRapMedHigh , 21, 1, 1);
book(_IncRapMedLow , 22, 1, 1);
book(_IncRapLow , 23, 1, 1);
}
void analyze(const Event& e) {
-
- // Get event weight for histo filling
- const double weight = 1.0;
-
-
// Final state of unstable particles to get particle spectra
const UnstableFinalState& ufs = apply<UnstableFinalState>(e, "UFS");
foreach (const Particle& p, ufs.particles()) {
if (p.abspid() != 443) continue;
const GenVertex* gv = p.genParticle()->production_vertex();
bool nonPrompt = false;
if (gv) {
foreach (const GenParticle* pi, Rivet::particles(gv, HepMC::ancestors)) {
const PdgId pid2 = pi->pdg_id();
if (PID::isHadron(pid2) && PID::hasBottom(pid2)) {
nonPrompt = true;
break;
}
}
}
double absrap = p.absrap();
double xp = p.perp();
if (absrap<=2.4 and absrap>2.) {
- if (nonPrompt) _nonPrRapHigh->fill(xp, weight);
- else if (!nonPrompt) _PrRapHigh->fill(xp, weight);
- _IncRapHigh->fill(xp, weight);
+ if (nonPrompt) _nonPrRapHigh->fill(xp);
+ else if (!nonPrompt) _PrRapHigh->fill(xp);
+ _IncRapHigh->fill(xp);
}
else if (absrap<=2. and absrap>1.5) {
- if (nonPrompt) _nonPrRapMedHigh->fill(xp, weight);
- else if (!nonPrompt) _PrRapMedHigh->fill(xp, weight);
- _IncRapMedHigh->fill(xp, weight);
+ if (nonPrompt) _nonPrRapMedHigh->fill(xp);
+ else if (!nonPrompt) _PrRapMedHigh->fill(xp);
+ _IncRapMedHigh->fill(xp);
}
else if (absrap<=1.5 and absrap>0.75) {
- if (nonPrompt) _nonPrRapMedLow->fill(xp, weight);
- else if (!nonPrompt) _PrRapMedLow->fill(xp, weight);
- _IncRapMedLow->fill(xp, weight);
+ if (nonPrompt) _nonPrRapMedLow->fill(xp);
+ else if (!nonPrompt) _PrRapMedLow->fill(xp);
+ _IncRapMedLow->fill(xp);
}
else if (absrap<=0.75) {
- if (nonPrompt) _nonPrRapLow->fill(xp, weight);
- else if (!nonPrompt) _PrRapLow->fill(xp, weight);
- _IncRapLow->fill(xp, weight);
+ if (nonPrompt) _nonPrRapLow->fill(xp);
+ else if (!nonPrompt) _PrRapLow->fill(xp);
+ _IncRapLow->fill(xp);
}
}
}
/// Finalize
void finalize() {
double factor = crossSection()/nanobarn*0.0593;
scale(_PrRapHigh , factor/sumOfWeights());
scale(_PrRapMedHigh , factor/sumOfWeights());
scale(_PrRapMedLow , factor/sumOfWeights());
scale(_PrRapLow , factor/sumOfWeights());
scale(_nonPrRapHigh , factor/sumOfWeights());
scale(_nonPrRapMedHigh, factor/sumOfWeights());
scale(_nonPrRapMedLow , factor/sumOfWeights());
scale(_nonPrRapLow , factor/sumOfWeights());
scale(_IncRapHigh , 1000.*factor/sumOfWeights());
scale(_IncRapMedHigh , 1000.*factor/sumOfWeights());
scale(_IncRapMedLow , 1000.*factor/sumOfWeights());
scale(_IncRapLow , 1000.*factor/sumOfWeights());
}
//@}
private:
Histo1DPtr _nonPrRapHigh;
Histo1DPtr _nonPrRapMedHigh;
Histo1DPtr _nonPrRapMedLow;
Histo1DPtr _nonPrRapLow;
Histo1DPtr _PrRapHigh;
Histo1DPtr _PrRapMedHigh;
Histo1DPtr _PrRapMedLow;
Histo1DPtr _PrRapLow;
Histo1DPtr _IncRapHigh;
Histo1DPtr _IncRapMedHigh;
Histo1DPtr _IncRapMedLow;
Histo1DPtr _IncRapLow;
//@}
};
// The hook for the plugin system
DECLARE_RIVET_PLUGIN(ATLAS_2011_S9035664);
}

File Metadata

Mime Type
text/x-diff
Expires
Tue, Nov 19, 6:29 PM (1 d, 21 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
3805596
Default Alt Text
(27 KB)

Event Timeline