Page MenuHomeHEPForge
Files F19250704

No OneTemporary
Actions

Size
26 KB
Referenced Files
None
Subscribers
None

View Options

diff --git a/src/Analyses/ATLAS_2011_I945498.cc b/src/Analyses/ATLAS_2011_I945498.cc
--- a/src/Analyses/ATLAS_2011_I945498.cc
+++ b/src/Analyses/ATLAS_2011_I945498.cc
@@ -1,302 +1,292 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/RivetYODA.hh"
#include "Rivet/Tools/Logging.hh"
#include "Rivet/Projections/ZFinder.hh"
#include "Rivet/Particle.fhh"
#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"
#include "Rivet/Projections/ClusteredPhotons.hh"
namespace Rivet {
class ATLAS_2011_I945498 : public Analysis {
public:
/// @name Constructors etc.
//@{
/// Constructor
ATLAS_2011_I945498()
: Analysis("ATLAS_2011_I945498")
{
- /// @todo Set whether your finalize method needs the generator cross section
setNeedsCrossSection(true);
for (size_t chn = 0; chn < 3; ++chn) {
weights_nj0[chn] = 0.0;
weights_nj1[chn] = 0.0;
weights_nj2[chn] = 0.0;
weights_nj3[chn] = 0.0;
weights_nj4[chn] = 0.0;
}
}
//@}
public:
-
/// Book histograms and initialise projections before the run
void init() {
// Set up projections
ZFinder zfinder_mu(-2.4, 2.4, 20, PID::MUON, 66.0*GeV, 116.0*GeV, 0.1, true, false);
addProjection(zfinder_mu, "ZFinder_mu");
std::vector<std::pair<double, double> > eta_e;
- eta_e.push_back(make_pair(-2.47,-1.52));
- eta_e.push_back(make_pair(-1.37,1.37));
- eta_e.push_back(make_pair(1.52,2.47));
+ eta_e.push_back(make_pair(-2.47, -1.52));
+ eta_e.push_back(make_pair(-1.37, 1.37));
+ eta_e.push_back(make_pair(1.52, 2.47));
ZFinder zfinder_el(eta_e, 20, PID::ELECTRON, 66.0*GeV, 116.0*GeV, 0.1, true, false);
addProjection(zfinder_el, "ZFinder_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);
FastJets jets(remfs, FastJets::ANTIKT, 0.4);
jets.useInvisibles();
addProjection(jets, "jets");
// 0=el, 1=mu, 2=comb
for (size_t chn = 0; chn < 3; ++chn) {
_h_njet_incl[chn] = bookHisto1D(1, 1, chn+1);
_h_njet_ratio[chn] = bookScatter2D(2, 1, chn+1);
_h_ptjet[chn] = bookHisto1D(3, 1, chn+1);
_h_ptlead[chn] = bookHisto1D(4, 1, chn+1);
_h_ptseclead[chn] = bookHisto1D(5, 1, chn+1);
_h_yjet[chn] = bookHisto1D(6, 1, chn+1);
_h_ylead[chn] = bookHisto1D(7, 1, chn+1);
_h_yseclead[chn] = bookHisto1D(8, 1, chn+1);
_h_mass[chn] = bookHisto1D(9, 1, chn+1);
_h_deltay[chn] = bookHisto1D(10, 1, chn+1);
_h_deltaphi[chn] = bookHisto1D(11, 1, chn+1);
_h_deltaR[chn] = bookHisto1D(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 Event& event) {
- FourMomentum l1=zf->constituents()[0].momentum();
- FourMomentum l2=zf->constituents()[1].momentum();
+ FourMomentum l1 = zf->constituents()[0].momentum();
+ FourMomentum l2 = zf->constituents()[1].momentum();
Jets jets;
-
foreach (const Jet& jet, applyProjection<FastJets>(event, "jets").jetsByPt(30.0*GeV)) {
FourMomentum jmom = jet.momentum();
if (fabs(jmom.rapidity()) < 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 = event.weight();
vector<const ZFinder*> zfs;
zfs.push_back(& (applyProjection<ZFinder>(event, "ZFinder_el")));
zfs.push_back(& (applyProjection<ZFinder>(event, "ZFinder_mu")));
+ // Require exactly one electronic or muonic Z-decay in the event
+ if (!( (zfs[0]->bosons().size() == 1 && zfs[1]->bosons().size() != 1) ||
+ (zfs[1]->bosons().size() == 1 && zfs[0]->bosons().size() != 1) )) vetoEvent;
+ int chn = zfs[0]->bosons().size() == 1 ? 0 : 1;
- // Require exactly one electronic or muonic Z-decay in the event
- if (!( (zfs[0]->bosons().size()==1 && zfs[1]->bosons().size()!=1) || (zfs[1]->bosons().size()==1 && zfs[0]->bosons().size()!=1) )) vetoEvent;
-
- int chn = (zfs[0]->bosons().size()==1 ? 0:1);
-
- Jets jets=selectJets(zfs[chn], event);
-
+ Jets jets = selectJets(zfs[chn], event);
/// TODO: Holger wrote in his commit message that the njet_ratio
/// histograms need fixing/checking, and therefore the analysis
/// is marked unvalidated!
// Some silly weight counters for the njet-ratio histo
// --- not sure about the njet=0 case, the Figure ca[ption says
// that selected events require at least one jet with 20 GeV
- if (jets.size() == 0) {
+ switch (jets.size()) {
+ case 0:
weights_nj0[chn] += 1.0;
- }
- else if (jets.size() == 1) {
+ break;
+ case 1:
weights_nj0[chn] += 1.0;
weights_nj1[chn] += 1.0;
- }
- else if (jets.size() == 2) {
+ break;
+ case 2:
weights_nj0[chn] += 1.0;
weights_nj1[chn] += 1.0;
weights_nj2[chn] += 1.0;
- }
- else if (jets.size() == 3) {
+ break;
+ case 3:
weights_nj0[chn] += 1.0;
weights_nj1[chn] += 1.0;
weights_nj2[chn] += 1.0;
weights_nj3[chn] += 1.0;
- }
- else if (jets.size() >= 4) {
+ break;
+ default: // >= 4
weights_nj0[chn] += 1.0;
weights_nj1[chn] += 1.0;
weights_nj2[chn] += 1.0;
weights_nj3[chn] += 1.0;
weights_nj4[chn] += 1.0;
}
- if (jets.size() <1) vetoEvent;
+ // Require at least one jet
+ if (jets.size() < 1) vetoEvent;
- _h_njet_incl[chn] ->fill(jets.size(), weight);
- _h_njet_incl[2] ->fill(jets.size(), weight);
+ // Fill jet multiplicities
+ _h_njet_incl[chn]->fill(jets.size(), weight);
+ _h_njet_incl[2]->fill(jets.size(), weight);
// Loop over selected jets, fill inclusive jet distributions
- for (unsigned int ijet=0; ijet<jets.size(); ijet++){
- _h_ptjet[chn] ->fill(jets[ijet].momentum().pT()/GeV, weight);
- _h_ptjet[2] ->fill(jets[ijet].momentum().pT()/GeV, weight);
- _h_yjet[chn] ->fill(fabs(jets[ijet].momentum().rapidity()), weight);
- _h_yjet[2] ->fill(fabs(jets[ijet].momentum().rapidity()), weight);
+ for (size_t ijet = 0; ijet < jets.size(); ++ijet) {
+ _h_ptjet[chn]->fill(jets[ijet].momentum().pT()/GeV, weight);
+ _h_ptjet[2] ->fill(jets[ijet].momentum().pT()/GeV, weight);
+ _h_yjet[chn] ->fill(fabs(jets[ijet].momentum().rapidity()), weight);
+ _h_yjet[2] ->fill(fabs(jets[ijet].momentum().rapidity()), weight);
}
- // The leading jet histos
- if (jets.size()>=1) {
- double ptlead = jets[0].momentum().pT()/GeV;
- double yabslead = fabs(jets[0].momentum().rapidity());
- _h_ptlead[chn] ->fill(ptlead, weight);
- _h_ptlead[2] ->fill(ptlead, weight);
+ // Leading jet histos
+ const double ptlead = jets[0].momentum().pT()/GeV;
+ const double yabslead = fabs(jets[0].momentum().rapidity());
+ _h_ptlead[chn]->fill(ptlead, weight);
+ _h_ptlead[2] ->fill(ptlead, weight);
+ _h_ylead[chn] ->fill(yabslead, weight);
+ _h_ylead[2] ->fill(yabslead, weight);
- _h_ylead[chn] ->fill(yabslead, weight);
- _h_ylead[2] ->fill(yabslead, weight);
+ if (jets.size() >= 2) {
+ // Second jet histos
+ const double pt2ndlead = jets[1].momentum().pT()/GeV;
+ const double yabs2ndlead = fabs(jets[1].momentum().rapidity());
+ _h_ptseclead[chn] ->fill(pt2ndlead, weight);
+ _h_ptseclead[2] ->fill(pt2ndlead, weight);
+ _h_yseclead[chn] ->fill(yabs2ndlead, weight);
+ _h_yseclead[2] ->fill(yabs2ndlead, weight);
+
+ // Dijet histos
+ const double deltaphi = fabs(deltaPhi(jets[1], jets[0]));
+ const double deltarap = fabs(jets[0].momentum().rapidity() - jets[1].momentum().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, weight);
+ _h_mass[2] ->fill(mass, weight);
+ _h_deltay[chn] ->fill(deltarap, weight);
+ _h_deltay[2] ->fill(deltarap, weight);
+ _h_deltaphi[chn]->fill(deltaphi, weight);
+ _h_deltaphi[2] ->fill(deltaphi, weight);
+ _h_deltaR[chn] ->fill(deltar, weight);
+ _h_deltaR[2] ->fill(deltar, weight);
}
- if (jets.size()>=2) {
- // The second to leading jet histos
- double pt2ndlead = jets[1].momentum().pT()/GeV;
- double yabs2ndlead = fabs(jets[1].momentum().rapidity());
-
- _h_ptseclead[chn] ->fill(pt2ndlead, weight);
- _h_ptseclead[2] ->fill(pt2ndlead, weight);
- _h_yseclead[chn] ->fill(yabs2ndlead, weight);
- _h_yseclead[2] ->fill(yabs2ndlead, weight);
-
- // Dijet histos
- double deltaphi = fabs(deltaPhi(jets[1], jets[0]));
- double deltarap = fabs(jets[0].momentum().rapidity() - jets[1].momentum().rapidity()) ;
- double deltar = fabs(deltaR(jets[0], jets[1], RAPIDITY));
- double mass = (jets[0].momentum() + jets[1].momentum()).mass();
-
- _h_mass[chn] ->fill(mass, weight);
- _h_mass[2] ->fill(mass, weight);
- _h_deltay[chn] ->fill(deltarap, weight);
- _h_deltay[2] ->fill(deltarap, weight);
- _h_deltaphi[chn] ->fill(deltaphi, weight);
- _h_deltaphi[2] ->fill(deltaphi, weight);
- _h_deltaR[chn] ->fill(deltar, weight);
- _h_deltaR[2] ->fill(deltar, weight);
- }
}
- /// Normalise histograms etc., after the run
+ /// @name Ratio calculator util functions
+ //@{
- std::vector<double> ratio(double a, double b) {
- double ratio=0.0;
- double ratio_err=0.0;
- std::vector<double> temp;
- cout << "a: " << a << " b: " << b << endl;
- if (b>0. && a>0.) {
- ratio = a/b;
- ratio_err = sqrt(a/pow(b,2) + pow(a,2)*b/pow(b,4));
- }
- temp.push_back(ratio);
- temp.push_back(ratio_err);
- return temp;
+ /// 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;
+ }
+
+ /// Calculate combined ratio from muon and electron channels
+ double comb_ratio(double* as, double* bs) {
+ return ratio(as[0]+as[1], bs[0]+bs[1]);
+ }
+
+ /// Calculate combined ratio error from muon and electron channels
+ double comb_ratio_err(double* as, double* bs) {
+ return ratio_err(as[0]+as[1], bs[0]+bs[1]);
+ }
+
+ //@}
+
+
void finalize() {
- // Fill RATIO histograms (DataPointSets)
- vector<double> yvals_ratio[3];
- vector<double> yerrs_ratio[3];
+ // Fill ratio histograms
+ for (size_t chn = 0; chn < 2; ++chn) {
+ _h_njet_ratio[chn]->point(0).setY( ratio(weights_nj1[chn], weights_nj0[chn]), ratio_err(weights_nj1[chn], weights_nj0[chn]) );
+ _h_njet_ratio[chn]->point(1).setY( ratio(weights_nj2[chn], weights_nj1[chn]), ratio_err(weights_nj2[chn], weights_nj1[chn]) );
+ _h_njet_ratio[chn]->point(2).setY( ratio(weights_nj3[chn], weights_nj2[chn]), ratio_err(weights_nj3[chn], weights_nj2[chn]) );
+ _h_njet_ratio[chn]->point(3).setY( ratio(weights_nj4[chn], weights_nj3[chn]), ratio_err(weights_nj4[chn], weights_nj3[chn]) );
+ }
+ _h_njet_ratio[2]->point(0).setY( comb_ratio(weights_nj1, weights_nj0), comb_ratio_err(weights_nj1, weights_nj0) );
+ _h_njet_ratio[2]->point(1).setY( comb_ratio(weights_nj2, weights_nj1), comb_ratio_err(weights_nj2, weights_nj1) );
+ _h_njet_ratio[2]->point(2).setY( comb_ratio(weights_nj3, weights_nj2), comb_ratio_err(weights_nj3, weights_nj2) );
+ _h_njet_ratio[2]->point(3).setY( comb_ratio(weights_nj4, weights_nj3), comb_ratio_err(weights_nj4, weights_nj3) );
- for (size_t chn = 0; chn < 2; ++chn) {
- yvals_ratio[chn].push_back(ratio(weights_nj1[chn], weights_nj0[chn])[0]);
- yvals_ratio[chn].push_back(ratio(weights_nj2[chn], weights_nj1[chn])[0]);
- yvals_ratio[chn].push_back(ratio(weights_nj3[chn], weights_nj2[chn])[0]);
- yvals_ratio[chn].push_back(ratio(weights_nj4[chn], weights_nj3[chn])[0]);
- // Errors
- yerrs_ratio[chn].push_back(ratio(weights_nj1[chn], weights_nj0[chn])[1]);
- yerrs_ratio[chn].push_back(ratio(weights_nj2[chn], weights_nj1[chn])[1]);
- yerrs_ratio[chn].push_back(ratio(weights_nj3[chn], weights_nj2[chn])[1]);
- yerrs_ratio[chn].push_back(ratio(weights_nj4[chn], weights_nj3[chn])[1]);
- // Actually fill histo
- /// \todo YODA
- //_h_njet_ratio[chn] ->setCoordinate(1, yvals_ratio[chn], yerrs_ratio[chn]);
- // Combined histos
- yvals_ratio[2].push_back(ratio(weights_nj1[chn], weights_nj0[chn])[0]);
- yvals_ratio[2].push_back(ratio(weights_nj2[chn], weights_nj1[chn])[0]);
- yvals_ratio[2].push_back(ratio(weights_nj3[chn], weights_nj2[chn])[0]);
- yvals_ratio[2].push_back(ratio(weights_nj4[chn], weights_nj3[chn])[0]);
- // Errors
- yerrs_ratio[2].push_back(ratio(weights_nj1[chn], weights_nj0[chn])[1]);
- yerrs_ratio[2].push_back(ratio(weights_nj2[chn], weights_nj1[chn])[1]);
- yerrs_ratio[2].push_back(ratio(weights_nj3[chn], weights_nj2[chn])[1]);
- yerrs_ratio[2].push_back(ratio(weights_nj4[chn], weights_nj3[chn])[1]);
- }
- /// \todo YODA
- //_h_njet_ratio[2] ->setCoordinate(1, yvals_ratio[2], yerrs_ratio[2]);
-
-
+ // Scale other histos
const double xs = crossSectionPerEvent()/picobarn;
for (size_t chn = 0; chn < 3; ++chn) {
scale(_h_njet_incl[chn], 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:
double weights_nj0[3];
double weights_nj1[3];
double weights_nj2[3];
double weights_nj3[3];
double 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];
};
- // This global object acts as a hook for the plugin system
+
DECLARE_RIVET_PLUGIN(ATLAS_2011_I945498);
+
+
}
diff --git a/src/Analyses/ATLAS_2012_I1083318.cc b/src/Analyses/ATLAS_2012_I1083318.cc
--- a/src/Analyses/ATLAS_2012_I1083318.cc
+++ b/src/Analyses/ATLAS_2012_I1083318.cc
@@ -1,270 +1,264 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/RivetYODA.hh"
#include "Rivet/Tools/Logging.hh"
#include "Rivet/Projections/IdentifiedFinalState.hh"
#include "Rivet/Projections/VetoedFinalState.hh"
#include "Rivet/Projections/MissingMomentum.hh"
#include "Rivet/Projections/FastJets.hh"
#include "Rivet/Projections/LeptonClusters.hh"
#include "Rivet/Projections/LeadingParticlesFinalState.hh"
namespace Rivet {
class ATLAS_2012_I1083318 : public Analysis {
public:
/// @name Constructors etc.
//@{
/// Constructor
ATLAS_2012_I1083318()
: Analysis("ATLAS_2012_I1083318")
{ }
//@}
public:
/// @name Analysis methods
//@{
/// Book histograms and initialise projections before the run
void init() {
FinalState fs;
IdentifiedFinalState allleptons;
allleptons.acceptIdPair(PID::ELECTRON);
allleptons.acceptIdPair(PID::MUON);
- std::vector<std::pair<double, double> > etaRanges;
+ vector< pair<double, double> > etaRanges;
etaRanges.push_back(make_pair(-2.5, 2.5));
LeptonClusters leptons(fs, allleptons,
0.1, true,
etaRanges, 20.0*GeV);
addProjection(leptons, "leptons");
// Leading neutrinos for Etmiss
LeadingParticlesFinalState neutrinos(fs);
neutrinos.addParticleIdPair(PID::NU_E);
neutrinos.addParticleIdPair(PID::NU_MU);
neutrinos.setLeadingOnly(true);
addProjection(neutrinos, "neutrinos");
// Input for the jets: "Neutrinos, electrons, and muons from decays of the
// massive W boson were not used"
VetoedFinalState veto;
veto.addVetoOnThisFinalState(leptons);
veto.addVetoOnThisFinalState(neutrinos);
FastJets jets(veto, FastJets::ANTIKT, 0.4);
jets.useInvisibles(true);
addProjection(jets, "jets");
- for (size_t i=0; i<2; ++i) {
+ for (size_t i = 0; i < 2; ++i) {
_h_NjetIncl[i] = bookHisto1D(1, 1, i+1);
_h_RatioNjetIncl[i] = bookScatter2D(2, 1, i+1);
_h_FirstJetPt_1jet[i] = bookHisto1D(3, 1, i+1);
_h_FirstJetPt_2jet[i] = bookHisto1D(4, 1, i+1);
_h_FirstJetPt_3jet[i] = bookHisto1D(5, 1, i+1);
_h_FirstJetPt_4jet[i] = bookHisto1D(6, 1, i+1);
_h_SecondJetPt_2jet[i] = bookHisto1D(7, 1, i+1);
_h_SecondJetPt_3jet[i] = bookHisto1D(8, 1, i+1);
_h_SecondJetPt_4jet[i] = bookHisto1D(9, 1, i+1);
_h_ThirdJetPt_3jet[i] = bookHisto1D(10, 1, i+1);
_h_ThirdJetPt_4jet[i] = bookHisto1D(11, 1, i+1);
_h_FourthJetPt_4jet[i] = bookHisto1D(12, 1, i+1);
_h_Ht_1jet[i] = bookHisto1D(13, 1, i+1);
_h_Ht_2jet[i] = bookHisto1D(14, 1, i+1);
_h_Ht_3jet[i] = bookHisto1D(15, 1, i+1);
_h_Ht_4jet[i] = bookHisto1D(16, 1, i+1);
_h_Minv_2jet[i] = bookHisto1D(17, 1, i+1);
_h_Minv_3jet[i] = bookHisto1D(18, 1, i+1);
_h_Minv_4jet[i] = bookHisto1D(19, 1, i+1);
_h_JetRapidity[i] = bookHisto1D(20, 1, i+1);
_h_DeltaYElecJet[i] = bookHisto1D(21, 1, i+1);
_h_SumYElecJet[i] = bookHisto1D(22, 1, i+1);
_h_DeltaR_2jet[i] = bookHisto1D(23, 1, i+1);
_h_DeltaY_2jet[i] = bookHisto1D(24, 1, i+1);
_h_DeltaPhi_2jet[i] = bookHisto1D(25, 1, i+1);
}
}
/// Perform the per-event analysis
void analyze(const Event& event) {
const double weight = event.weight();
const vector<ClusteredLepton>& leptons = applyProjection<LeptonClusters>(event, "leptons").clusteredLeptons();
Particles neutrinos = applyProjection<FinalState>(event, "neutrinos").particlesByPt();
- if (leptons.size()!=1 || (neutrinos.size()==0)) {
+ if (leptons.size() != 1 || (neutrinos.size() == 0)) {
vetoEvent;
}
- FourMomentum lepton=leptons[0].momentum();
+ FourMomentum lepton = leptons[0].momentum();
FourMomentum p_miss = neutrinos[0].momentum();
- if (p_miss.Et()<25.0*GeV) {
+ if (p_miss.Et() < 25.0*GeV) {
vetoEvent;
}
- double mT=sqrt(2.0*lepton.pT()*p_miss.Et()*(1.0-cos(lepton.phi()-p_miss.phi())));
- if (mT<40.0*GeV) {
+ double mT = sqrt(2.0 * lepton.pT() * p_miss.Et() * (1.0 - cos( lepton.phi()-p_miss.phi()) ) );
+ if (mT < 40.0*GeV) {
vetoEvent;
}
- double jetcuts[] = {30.0*GeV, 20.0*GeV};
+ double jetcuts[] = { 30.0*GeV, 20.0*GeV };
const FastJets& jetpro = applyProjection<FastJets>(event, "jets");
- for (size_t i=0; i<2; ++i) {
+ for (size_t i = 0; i < 2; ++i) {
vector<FourMomentum> jets;
- double HT=lepton.pT()+p_miss.pT();
+ double HT = lepton.pT() + p_miss.pT();
foreach (const Jet& jet, jetpro.jetsByPt(jetcuts[i])) {
- if (fabs(jet.momentum().rapidity())<4.4 &&
- deltaR(lepton, jet.momentum())>0.5) {
+ if (fabs(jet.momentum().rapidity()) < 4.4 && deltaR(lepton, jet.momentum()) > 0.5) {
jets.push_back(jet.momentum());
HT += jet.momentum().pT();
}
}
_h_NjetIncl[i]->fill(0.0, weight);
// Njet>=1 observables
- if (jets.size()<1) continue;
+ if (jets.size() < 1) continue;
_h_NjetIncl[i]->fill(1.0, weight);
_h_FirstJetPt_1jet[i]->fill(jets[0].pT(), weight);
_h_JetRapidity[i]->fill(jets[0].rapidity(), weight);
_h_Ht_1jet[i]->fill(HT, weight);
_h_DeltaYElecJet[i]->fill(lepton.rapidity()-jets[0].rapidity(), weight);
_h_SumYElecJet[i]->fill(lepton.rapidity()+jets[0].rapidity(), weight);
// Njet>=2 observables
- if (jets.size()<2) continue;
+ if (jets.size() < 2) continue;
_h_NjetIncl[i]->fill(2.0, weight);
_h_FirstJetPt_2jet[i]->fill(jets[0].pT(), weight);
_h_SecondJetPt_2jet[i]->fill(jets[1].pT(), weight);
_h_Ht_2jet[i]->fill(HT, weight);
double m2_2jet = FourMomentum(jets[0]+jets[1]).mass2();
_h_Minv_2jet[i]->fill(m2_2jet>0.0 ? sqrt(m2_2jet) : 0.0, weight);
_h_DeltaR_2jet[i]->fill(deltaR(jets[0], jets[1]), weight);
_h_DeltaY_2jet[i]->fill(jets[0].rapidity()-jets[1].rapidity(), weight);
_h_DeltaPhi_2jet[i]->fill(deltaPhi(jets[0], jets[1]), weight);
// Njet>=3 observables
- if (jets.size()<3) continue;
+ if (jets.size() < 3) continue;
_h_NjetIncl[i]->fill(3.0, weight);
_h_FirstJetPt_3jet[i]->fill(jets[0].pT(), weight);
_h_SecondJetPt_3jet[i]->fill(jets[1].pT(), weight);
_h_ThirdJetPt_3jet[i]->fill(jets[2].pT(), weight);
_h_Ht_3jet[i]->fill(HT, weight);
double m2_3jet = FourMomentum(jets[0]+jets[1]+jets[2]).mass2();
_h_Minv_3jet[i]->fill(m2_3jet>0.0 ? sqrt(m2_3jet) : 0.0, weight);
// Njet>=4 observables
- if (jets.size()<4) continue;
+ if (jets.size() < 4) continue;
_h_NjetIncl[i]->fill(4.0, weight);
_h_FirstJetPt_4jet[i]->fill(jets[0].pT(), weight);
_h_SecondJetPt_4jet[i]->fill(jets[1].pT(), weight);
_h_ThirdJetPt_4jet[i]->fill(jets[2].pT(), weight);
_h_FourthJetPt_4jet[i]->fill(jets[3].pT(), weight);
_h_Ht_4jet[i]->fill(HT, weight);
double m2_4jet = FourMomentum(jets[0]+jets[1]+jets[2]+jets[3]).mass2();
_h_Minv_4jet[i]->fill(m2_4jet>0.0 ? sqrt(m2_4jet) : 0.0, weight);
// Njet>=5 observables
- if (jets.size()<5) continue;
+ if (jets.size() < 5) continue;
_h_NjetIncl[i]->fill(5.0, weight);
}
}
/// Normalise histograms etc., after the run
void finalize() {
- for (size_t i=0; i<2; ++i) {
- /// @todo YODA
- //// first construct jet multi ratio
- //int Nbins = _h_NjetIncl[i]->numBins();
- //std::vector<double> ratio(Nbins-1, 0.0);
- //std::vector<double> err(Nbins-1, 0.0);
- //for (int n = 0; n < Nbins-1; ++n) {
- // if (_h_NjetIncl[i]->binHeight(n) > 0.0 && _h_NjetIncl[i]->binHeight(n+1) > 0.0) {
- // ratio[n] = _h_NjetIncl[i]->binHeight(n+1)/_h_NjetIncl[i]->binHeight(n);
- // double relerr_n = _h_NjetIncl[i]->binError(n)/_h_NjetIncl[i]->binHeight(n);
- // double relerr_m = _h_NjetIncl[i]->binError(n+1)/_h_NjetIncl[i]->binHeight(n+1);
- // err[n] = ratio[n] * (relerr_n + relerr_m);
- // }
- //}
- //_h_RatioNjetIncl[i]->setCoordinate(1, ratio, err);
+ for (size_t i = 0; i < 2; ++i) {
- // scale all histos to the cross section
- double factor = crossSection()/sumOfWeights();
+ // Construct jet multiplicity ratio
+ for (int n = 0; n < _h_NjetIncl[i]->numPoints()-1; ++n) {
+ YODA::Bin& b0 = _h_NjetIncl[i]->bin(n);
+ YODA::Bin& b1 = _h_NjetIncl[i]->bin(n+1);
+ if (b0.height() == 0.0 || b1.height() == 0.0) continue;
+ _h_RatioNjetIncl[i]->point(n).setY(b1.height()/b0.height());
+ _h_RatioNjetIncl[i]->point(n).setYErr(b1.height()/b0.height() * (b0.relErr() + b1.relErr()));
+ }
+
+ // Scale all histos to the cross section
+ const double factor = crossSection()/sumOfWeights();
scale(_h_DeltaPhi_2jet[i], factor);
scale(_h_DeltaR_2jet[i], factor);
scale(_h_DeltaY_2jet[i], factor);
scale(_h_DeltaYElecJet[i], factor);
scale(_h_FirstJetPt_1jet[i], factor);
scale(_h_FirstJetPt_2jet[i], factor);
scale(_h_FirstJetPt_3jet[i], factor);
scale(_h_FirstJetPt_4jet[i], factor);
scale(_h_FourthJetPt_4jet[i], factor);
scale(_h_Ht_1jet[i], factor);
scale(_h_Ht_2jet[i], factor);
scale(_h_Ht_3jet[i], factor);
scale(_h_Ht_4jet[i], factor);
scale(_h_JetRapidity[i], factor);
scale(_h_Minv_2jet[i], factor);
scale(_h_Minv_3jet[i], factor);
scale(_h_Minv_4jet[i], factor);
scale(_h_NjetIncl[i], factor);
scale(_h_SecondJetPt_2jet[i], factor);
scale(_h_SecondJetPt_3jet[i], factor);
scale(_h_SecondJetPt_4jet[i], factor);
scale(_h_SumYElecJet[i], factor);
scale(_h_ThirdJetPt_3jet[i], factor);
scale(_h_ThirdJetPt_4jet[i], factor);
}
}
//@}
private:
/// @name Histograms
//@{
Histo1DPtr _h_DeltaPhi_2jet[2];
Histo1DPtr _h_DeltaR_2jet[2];
Histo1DPtr _h_DeltaY_2jet[2];
Histo1DPtr _h_DeltaYElecJet[2];
Histo1DPtr _h_FirstJetPt_1jet[2];
Histo1DPtr _h_FirstJetPt_2jet[2];
Histo1DPtr _h_FirstJetPt_3jet[2];
Histo1DPtr _h_FirstJetPt_4jet[2];
Histo1DPtr _h_FourthJetPt_4jet[2];
Histo1DPtr _h_Ht_1jet[2];
Histo1DPtr _h_Ht_2jet[2];
Histo1DPtr _h_Ht_3jet[2];
Histo1DPtr _h_Ht_4jet[2];
Histo1DPtr _h_JetRapidity[2];
Histo1DPtr _h_Minv_2jet[2];
Histo1DPtr _h_Minv_3jet[2];
Histo1DPtr _h_Minv_4jet[2];
Histo1DPtr _h_NjetIncl[2];
Scatter2DPtr _h_RatioNjetIncl[2];
Histo1DPtr _h_SecondJetPt_2jet[2];
Histo1DPtr _h_SecondJetPt_3jet[2];
Histo1DPtr _h_SecondJetPt_4jet[2];
Histo1DPtr _h_SumYElecJet[2];
Histo1DPtr _h_ThirdJetPt_3jet[2];
Histo1DPtr _h_ThirdJetPt_4jet[2];
//@}
};
// The hook for the plugin system
DECLARE_RIVET_PLUGIN(ATLAS_2012_I1083318);
}

File Metadata

Mime Type
text/x-diff
Expires
Tue, Sep 30, 5:43 AM (5 h, 49 m)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
6566298
Default Alt Text
(26 KB)

Event Timeline

HEPForge.org · .