Page Menu
Home
HEPForge
Search
Configure Global Search
Log In
Files
F7878625
No One
Temporary
Actions
View File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Flag For Later
Size
27 KB
Subscribers
None
View Options
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
Details
Attached
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)
Attached To
rRIVETHG rivethg
Event Timeline
Log In to Comment