Page MenuHomeHEPForge
Files F7878989

No OneTemporary
Actions

View Options

diff --git a/analyses/pluginLHCb/LHCB_2010_I867355.cc b/analyses/pluginLHCb/LHCB_2010_I867355.cc
--- a/analyses/pluginLHCb/LHCB_2010_I867355.cc
+++ b/analyses/pluginLHCb/LHCB_2010_I867355.cc
@@ -1,90 +1,88 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Particle.hh"
namespace Rivet {
class LHCB_2010_I867355 : public Analysis {
public:
LHCB_2010_I867355() : Analysis("LHCB_2010_I867355")
{ }
void init() {
//@ Results are presented for two different fragmentation functions, LEP and Tevatron. Therefore, we have two sets of histograms.
book(_h_sigma_vs_eta_lep ,1, 1, 1);
book(_h_sigma_vs_eta_tvt ,1, 1, 2);
book(_h_sigma_total_lep ,2, 1, 1);
book(_h_sigma_total_tvt ,2, 1, 2);
}
/// Perform the per-event analysis
void analyze(const Event& event) {
- double weight = 1.0;
-
Particles bhadrons;
foreach (const GenParticle* p, particles(event.genEvent())) {
if (!( PID::isHadron( p->pdg_id() ) && PID::hasBottom( p->pdg_id() )) ) continue;
const GenVertex* dv = p->end_vertex();
bool hasBdaughter = false;
if ( PID::isHadron( p->pdg_id() ) && PID::hasBottom( p->pdg_id() )) { // selecting b-hadrons
if (dv) {
for (GenVertex::particles_out_const_iterator pp = dv->particles_out_const_begin() ; pp != dv->particles_out_const_end() ; ++pp) {
if (PID::isHadron( (*pp)->pdg_id() ) && PID::hasBottom( (*pp)->pdg_id() )) {
hasBdaughter = true;
}
}
}
}
if (hasBdaughter) continue; // continue if the daughter is another b-hadron
bhadrons += Particle(*p);
}
foreach (const Particle& particle, bhadrons) {
// take fabs() to use full statistics and then multiply weight by 0.5 because LHCb is single-sided
double eta = fabs(particle.eta());
- _h_sigma_vs_eta_lep->fill( eta, 0.5*weight );
- _h_sigma_vs_eta_tvt->fill( eta, 0.5*weight );
+ _h_sigma_vs_eta_lep->fill( eta, 0.5 );
+ _h_sigma_vs_eta_tvt->fill( eta, 0.5 );
- _h_sigma_total_lep->fill( eta, 0.5*weight ); // histogram for full kinematic range
- _h_sigma_total_tvt->fill( eta, 0.5*weight ); // histogram for full kinematic range
+ _h_sigma_total_lep->fill( eta, 0.5 ); // histogram for full kinematic range
+ _h_sigma_total_tvt->fill( eta, 0.5 ); // histogram for full kinematic range
}
}
void finalize() {
double norm = crossSection()/microbarn/sumOfWeights();
double binwidth = 4.; // integrated over full rapidity space from 2 to 6.
// to get the avergae of b and bbar, we scale with 0.5
scale(_h_sigma_vs_eta_lep, 0.5*norm);
scale(_h_sigma_vs_eta_tvt, 0.5*norm);
scale(_h_sigma_total_lep, 0.5*norm*binwidth);
scale(_h_sigma_total_tvt, 0.5*norm*binwidth);
}
private:
Histo1DPtr _h_sigma_total_lep;
Histo1DPtr _h_sigma_total_tvt;
Histo1DPtr _h_sigma_vs_eta_lep;
Histo1DPtr _h_sigma_vs_eta_tvt;
};
// Hook for the plugin system
DECLARE_RIVET_PLUGIN(LHCB_2010_I867355);
}
diff --git a/analyses/pluginLHCb/LHCB_2010_S8758301.cc b/analyses/pluginLHCb/LHCB_2010_S8758301.cc
--- a/analyses/pluginLHCb/LHCB_2010_S8758301.cc
+++ b/analyses/pluginLHCb/LHCB_2010_S8758301.cc
@@ -1,341 +1,554 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/UnstableFinalState.hh"
#include "Rivet/Math/Constants.hh"
#include "Rivet/Math/Units.hh"
#include "HepMC/GenEvent.h"
#include "HepMC/GenParticle.h"
#include "HepMC/GenVertex.h"
#include "HepMC/SimpleVector.h"
namespace Rivet {
using namespace HepMC;
-using namespace std;
// Lifetime cut: longest living ancestor ctau < 10^-11 [m]
namespace {
const double MAX_CTAU = 1.0E-11; // [m]
const double MIN_PT = 0.0001; // [GeV/c]
}
class LHCB_2010_S8758301 : public Analysis {
public:
/// @name Constructors etc.
//@{
/// Constructor
LHCB_2010_S8758301()
: Analysis("LHCB_2010_S8758301"),
- sumKs0_30(0.0), sumKs0_35(0.0),
- sumKs0_40(0.0), sumKs0_badnull(0),
+ sumKs0_badnull(0),
sumKs0_badlft(0), sumKs0_all(0),
sumKs0_outup(0), sumKs0_outdwn(0),
sum_low_pt_loss(0), sum_high_pt_loss(0)
{ }
//@}
/// @name Analysis methods
//@{
/// Book histograms and initialise projections before the run
void init() {
MSG_DEBUG("Initializing analysis!");
fillMap(partLftMap);
book(_h_K0s_pt_30 ,1,1,1);
book(_h_K0s_pt_35 ,1,1,2);
book(_h_K0s_pt_40 ,1,1,3);
book(_h_K0s_pt_y_30 ,2,1,1);
book(_h_K0s_pt_y_35 ,2,1,2);
book(_h_K0s_pt_y_40 ,2,1,3);
book(_h_K0s_pt_y_all ,3,1,1);
+
+ book(sumKs0_30, "TMP/sumKs0_30");
+ book(sumKs0_35, "TMP/sumKs0_35");
+ book(sumKs0_40, "TMP/sumKs0_40");
declare(UnstableFinalState(), "UFS");
}
/// Perform the per-event analysis
void analyze(const Event& event) {
int id;
double y, pT;
- const double weight = 1.0;
const UnstableFinalState& ufs = apply<UnstableFinalState>(event, "UFS");
double ancestor_lftime;
foreach (const Particle& p, ufs.particles()) {
id = p.pid();
if ((id != 310) && (id != -310)) continue;
sumKs0_all ++;
ancestor_lftime = 0.;
const GenParticle* long_ancestor = getLongestLivedAncestor(p, ancestor_lftime);
if ( !(long_ancestor) ) {
sumKs0_badnull ++;
continue;
}
if ( ancestor_lftime > MAX_CTAU ) {
sumKs0_badlft ++;
MSG_DEBUG("Ancestor " << long_ancestor->pdg_id() << ", ctau: " << ancestor_lftime << " [m]");
continue;
}
const FourMomentum& qmom = p.momentum();
y = 0.5 * log((qmom.E() + qmom.pz())/(qmom.E() - qmom.pz()));
pT = sqrt((qmom.px() * qmom.px()) + (qmom.py() * qmom.py()));
if (pT < MIN_PT) {
sum_low_pt_loss ++;
MSG_DEBUG("Small pT K^0_S: " << pT << " GeV/c.");
}
if (pT > 1.6) {
sum_high_pt_loss ++;
}
if (y > 2.5 && y < 4.0) {
- _h_K0s_pt_y_all->fill(pT, weight);
+ _h_K0s_pt_y_all->fill(pT);
if (y > 2.5 && y < 3.0) {
- _h_K0s_pt_y_30->fill(pT, weight);
- _h_K0s_pt_30->fill(pT, weight);
- sumKs0_30 += weight;
+ _h_K0s_pt_y_30->fill(pT);
+ _h_K0s_pt_30->fill(pT);
+ sumKs0_30->fill();
} else if (y > 3.0 && y < 3.5) {
- _h_K0s_pt_y_35->fill(pT, weight);
- _h_K0s_pt_35->fill(pT, weight);
- sumKs0_35 += weight;
+ _h_K0s_pt_y_35->fill(pT);
+ _h_K0s_pt_35->fill(pT);
+ sumKs0_35->fill();
} else if (y > 3.5 && y < 4.0) {
- _h_K0s_pt_y_40->fill(pT, weight);
- _h_K0s_pt_40->fill(pT, weight);
- sumKs0_40 += weight;
+ _h_K0s_pt_y_40->fill(pT);
+ _h_K0s_pt_40->fill(pT);
+ sumKs0_40->fill();
}
} else if (y < 2.5) {
sumKs0_outdwn ++;
} else if (y > 4.0) {
sumKs0_outup ++;
}
}
}
/// Normalise histograms etc., after the run
void finalize() {
MSG_DEBUG("Total number Ks0: " << sumKs0_all << endl
<< "Sum of weights: " << sumOfWeights() << endl
- << "Weight Ks0 (2.5 < y < 3.0): " << sumKs0_30 << endl
- << "Weight Ks0 (3.0 < y < 3.5): " << sumKs0_35 << endl
- << "Weight Ks0 (3.5 < y < 4.0): " << sumKs0_40 << endl
+ << "Weight Ks0 (2.5 < y < 3.0): " << double(sumKs0_30) << endl
+ << "Weight Ks0 (3.0 < y < 3.5): " << double(sumKs0_35) << endl
+ << "Weight Ks0 (3.5 < y < 4.0): " << double(sumKs0_40) << endl
<< "Nb. unprompt Ks0 [null mother]: " << sumKs0_badnull << endl
<< "Nb. unprompt Ks0 [mother lifetime exceeded]: " << sumKs0_badlft << endl
<< "Nb. Ks0 (y > 4.0): " << sumKs0_outup << endl
<< "Nb. Ks0 (y < 2.5): " << sumKs0_outdwn << endl
<< "Nb. Ks0 (pT < " << (MIN_PT/MeV) << " MeV/c): " << sum_low_pt_loss << endl
<< "Nb. Ks0 (pT > 1.6 GeV/c): " << sum_high_pt_loss << endl
<< "Cross-section [mb]: " << crossSection()/millibarn << endl
<< "Nb. events: " << numEvents());
// Compute cross-section; multiply by bin width for correct scaling
// cross-section given by Rivet in pb
double xsection_factor = crossSection()/sumOfWeights();
// Multiply bin width for correct scaling, xsection in mub
scale(_h_K0s_pt_30, 0.2*xsection_factor/microbarn);
scale(_h_K0s_pt_35, 0.2*xsection_factor/microbarn);
scale(_h_K0s_pt_40, 0.2*xsection_factor/microbarn);
// Divide by dy (rapidity window width), xsection in mb
scale(_h_K0s_pt_y_30, xsection_factor/0.5/millibarn);
scale(_h_K0s_pt_y_35, xsection_factor/0.5/millibarn);
scale(_h_K0s_pt_y_40, xsection_factor/0.5/millibarn);
scale(_h_K0s_pt_y_all, xsection_factor/1.5/millibarn);
}
//@}
private:
/// Get particle lifetime from hardcoded data
double getLifeTime(int pid) {
double lft = -1.0;
if (pid < 0) pid = - pid;
// Correct Pythia6 PIDs for f0(980), f0(1370) mesons
if (pid == 10331) pid = 30221;
if (pid == 10221) pid = 9010221;
map<int, double>::iterator pPartLft = partLftMap.find(pid);
// search stable particle list
if (pPartLft == partLftMap.end()) {
if (pid <= 100) return 0.0;
- for (unsigned int i=0; i < sizeof(stablePDGIds)/sizeof(unsigned int); i++ ) {
- if (pid == stablePDGIds[i]) { lft = 0.0; break; }
+ for ( auto id : stablePDGIds ) {
+ if (pid == id) { lft = 0.0; break; }
}
} else {
lft = (*pPartLft).second;
}
if (lft < 0.0)
MSG_ERROR("Could not determine lifetime for particle with PID " << pid
<< "... This K_s^0 will be considered unprompt!");
return lft;
}
const GenParticle* getLongestLivedAncestor(const Particle& p, double& lifeTime) {
- const GenParticle* ret = NULL;
+ const GenParticle* ret = nullptr;
lifeTime = 1.;
- if (p.genParticle() == NULL) return NULL;
+ if (p.genParticle() == nullptr) return nullptr;
const GenParticle* pmother = p.genParticle();
double longest_ctau = 0.;
double mother_ctau;
int mother_pid, n_inparts;
const GenVertex* ivertex = pmother->production_vertex();
while (ivertex) {
n_inparts = ivertex->particles_in_size();
- if (n_inparts < 1) {ret = NULL; break;} // error: should never happen!
- const GenVertex::particles_in_const_iterator iPart_invtx = ivertex->particles_in_const_begin();
+ if (n_inparts < 1) {ret = nullptr; break;} // error: should never happen!
+ const auto iPart_invtx = ivertex->particles_in_const_begin();
pmother = (*iPart_invtx); // first mother particle
mother_pid = pmother->pdg_id();
ivertex = pmother->production_vertex(); // get next vertex
if ( (mother_pid == 2212) || (mother_pid <= 100) ) {
- if (ret == NULL) ret = pmother;
+ if (ret == nullptr) ret = pmother;
continue;
}
mother_ctau = getLifeTime(mother_pid);
- if (mother_ctau < 0.) { ret= NULL; break; } // error:should never happen!
+ if (mother_ctau < 0.) { ret= nullptr; break; } // error:should never happen!
if (mother_ctau > longest_ctau) {
longest_ctau = mother_ctau;
ret = pmother;
}
}
if (ret) lifeTime = longest_ctau * c_light;
return ret;
}
// Fill the PDG Id to Lifetime[seconds] map
// Data was extract from LHCb Particle Table using ParticleSvc
bool fillMap(map<int, double> &m) {
- m[6] = 4.707703E-25; m[11] = 1.E+16; m[12] = 1.E+16;
- m[13] = 2.197019E-06; m[14] = 1.E+16; m[15] = 2.906E-13; m[16] = 1.E+16; m[22] = 1.E+16;
- m[23] = 2.637914E-25; m[24] = 3.075758E-25; m[25] = 9.4E-26; m[35] = 9.4E-26;
- m[36] = 9.4E-26; m[37] = 9.4E-26; m[84] = 3.335641E-13; m[85] = 1.290893E-12;
- m[111] = 8.4E-17; m[113] = 4.405704E-24; m[115] = 6.151516E-24; m[117] = 4.088275E-24;
- m[119] = 2.102914E-24; m[130] = 5.116E-08; m[150] = 1.525E-12; m[211] = 2.6033E-08;
- m[213] = 4.405704E-24; m[215] = 6.151516E-24; m[217] = 4.088275E-24; m[219] = 2.102914E-24;
- m[221] = 5.063171E-19; m[223] = 7.752794E-23; m[225] = 3.555982E-24; m[227] = 3.91793E-24;
- m[229] = 2.777267E-24; m[310] = 8.953E-11; m[313] = 1.308573E-23; m[315] = 6.038644E-24;
- m[317] = 4.139699E-24; m[319] = 3.324304E-24; m[321] = 1.238E-08; m[323] = 1.295693E-23;
- m[325] = 6.682357E-24; m[327] = 4.139699E-24; m[329] = 3.324304E-24; m[331] = 3.210791E-21;
- m[333] = 1.545099E-22; m[335] = 9.016605E-24; m[337] = 7.565657E-24; m[350] = 1.407125E-12;
- m[411] = 1.04E-12; m[413] = 6.856377E-21; m[415] = 1.778952E-23; m[421] = 4.101E-13;
- m[423] = 1.000003E-19; m[425] = 1.530726E-23; m[431] = 5.E-13; m[433] = 1.000003E-19;
- m[435] = 3.291061E-23; m[441] = 2.465214E-23; m[443] = 7.062363E-21; m[445] = 3.242425E-22;
- m[510] = 1.525E-12; m[511] = 1.525E-12; m[513] = 1.000019E-19; m[515] = 1.31E-23;
- m[521] = 1.638E-12; m[523] = 1.000019E-19; m[525] = 1.31E-23; m[530] = 1.536875E-12;
- m[531] = 1.472E-12; m[533] = 1.E-19; m[535] = 1.31E-23; m[541] = 4.5E-13;
- m[553] = 1.218911E-20; m[1112] = 4.539394E-24; m[1114] = 5.578069E-24; m[1116] = 1.994582E-24;
- m[1118] = 2.269697E-24; m[1212] = 4.539394E-24; m[1214] = 5.723584E-24; m[1216] = 1.994582E-24;
- m[1218] = 1.316424E-24; m[2112] = 8.857E+02; m[2114] = 5.578069E-24; m[2116] = 4.388081E-24;
- m[2118] = 2.269697E-24; m[2122] = 4.539394E-24; m[2124] = 5.723584E-24; m[2126] = 1.994582E-24;
- m[2128] = 1.316424E-24; m[2212] = 1.E+16; m[2214] = 5.578069E-24; m[2216] = 4.388081E-24;
- m[2218] = 2.269697E-24; m[2222] = 4.539394E-24; m[2224] = 5.578069E-24; m[2226] = 1.994582E-24;
- m[2228] = 2.269697E-24; m[3112] = 1.479E-10; m[3114] = 1.670589E-23; m[3116] = 5.485102E-24;
- m[3118] = 3.656734E-24; m[3122] = 2.631E-10; m[3124] = 4.219309E-23; m[3126] = 8.227653E-24;
- m[3128] = 3.291061E-24; m[3212] = 7.4E-20; m[3214] = 1.828367E-23; m[3216] = 5.485102E-24;
- m[3218] = 3.656734E-24; m[3222] = 8.018E-11; m[3224] = 1.838582E-23; m[3226] = 5.485102E-24;
- m[3228] = 3.656734E-24; m[3312] = 1.639E-10; m[3314] = 6.648608E-23; m[3322] = 2.9E-10;
- m[3324] = 7.233101E-23; m[3334] = 8.21E-11; m[4112] = 2.991874E-22; m[4114] = 4.088274E-23;
- m[4122] = 2.E-13; m[4132] = 1.12E-13; m[4212] = 3.999999E-22; m[4214] = 3.291061E-22;
- m[4222] = 2.951624E-22; m[4224] = 4.417531E-23; m[4232] = 4.42E-13; m[4332] = 6.9E-14;
- m[4412] = 3.335641E-13; m[4422] = 3.335641E-13; m[4432] = 3.335641E-13; m[5112] = 1.E-19;
- m[5122] = 1.38E-12; m[5132] = 1.42E-12; m[5142] = 1.290893E-12; m[5212] = 1.E-19;
- m[5222] = 1.E-19; m[5232] = 1.42E-12; m[5242] = 1.290893E-12; m[5312] = 1.E-19;
- m[5322] = 1.E-19; m[5332] = 1.55E-12; m[5342] = 1.290893E-12; m[5442] = 1.290893E-12;
- m[5512] = 1.290893E-12; m[5522] = 1.290893E-12; m[5532] = 1.290893E-12; m[5542] = 1.290893E-12;
- m[10111] = 2.48382E-24; m[10113] = 4.635297E-24; m[10115] = 2.54136E-24; m[10211] = 2.48382E-24;
- m[10213] = 4.635297E-24; m[10215] = 2.54136E-24; m[10223] = 1.828367E-24; m[10225] = 3.636531E-24;
- m[10311] = 2.437823E-24; m[10313] = 7.313469E-24; m[10315] = 3.538775E-24;
- m[10321] = 2.437823E-24; m[10323] = 7.313469E-24; m[10325] = 3.538775E-24;
- m[10331] = 4.804469E-24; m[10411] = 4.38E-24; m[10413] = 3.29E-23; m[10421] = 4.38E-24;
- m[10423] = 3.22653E-23; m[10431] = 6.5821E-22; m[10433] = 6.5821E-22; m[10441] = 6.453061E-23;
- m[10511] = 4.39E-24; m[10513] = 1.65E-23; m[10521] = 4.39E-24; m[10523] = 1.65E-23;
- m[10531] = 4.39E-24; m[10533] = 1.65E-23; m[11114] = 2.194041E-24; m[11116] = 1.828367E-24;
- m[11212] = 1.880606E-24; m[11216] = 1.828367E-24; m[12112] = 2.194041E-24;
- m[12114] = 2.194041E-24; m[12116] = 5.063171E-24; m[12126] = 1.828367E-24;
- m[12212] = 2.194041E-24; m[12214] = 2.194041E-24; m[12216] = 5.063171E-24;
- m[12224] = 2.194041E-24; m[12226] = 1.828367E-24; m[13112] = 6.582122E-24; m[13114] = 1.09702E-23;
- m[13116] = 5.485102E-24; m[13122] = 1.316424E-23; m[13124] = 1.09702E-23; m[13126] = 6.928549E-24;
- m[13212] = 6.582122E-24; m[13214] = 1.09702E-23; m[13216] = 5.485102E-24; m[13222] = 6.582122E-24;
- m[13224] = 1.09702E-23; m[13226] = 5.485102E-24;
- m[13312] = 4.135667E-22; m[13314] = 2.742551E-23;
- m[13324] = 2.742551E-23; m[14122] = 1.828367E-22; m[20022] = 1.E+16; m[20113] = 1.567172E-24;
- m[20213] = 1.567172E-24; m[20223] = 2.708692E-23; m[20313] = 3.782829E-24;
- m[20315] = 2.384827E-24; m[20323] = 3.782829E-24; m[20325] = 2.384827E-24;
- m[20333] = 1.198929E-23; m[20413] = 2.63E-24; m[20423] = 2.63E-24; m[20433] = 6.5821E-22;
- m[20443] = 7.395643E-22; m[20513] = 2.63E-24; m[20523] = 2.63E-24; m[20533] = 2.63E-24;
- m[21112] = 2.632849E-24; m[21114] = 3.291061E-24; m[21212] = 2.632849E-24;
- m[21214] = 6.582122E-24; m[22112] = 4.388081E-24; m[22114] = 3.291061E-24;
- m[22122] = 2.632849E-24; m[22124] = 6.582122E-24; m[22212] = 4.388081E-24;
- m[22214] = 3.291061E-24; m[22222] = 2.632849E-24; m[22224] = 3.291061E-24;
- m[23112] = 7.313469E-24; m[23114] = 2.991874E-24; m[23122] = 4.388081E-24;
- m[23124] = 6.582122E-24; m[23126] = 3.291061E-24; m[23212] = 7.313469E-24;
- m[23214] = 2.991874E-24; m[23222] = 7.313469E-24; m[23224] = 2.991874E-24;
- m[30113] = 2.632849E-24; m[30213] = 2.632849E-24; m[30221] = 1.880606E-24;
- m[30223] = 2.089563E-24; m[30313] = 2.056913E-24; m[30323] = 2.056913E-24;
- m[30443] = 2.419898E-23; m[31114] = 1.880606E-24; m[31214] = 3.291061E-24;
- m[32112] = 3.989164E-24; m[32114] = 1.880606E-24; m[32124] = 3.291061E-24;
- m[32212] = 3.989164E-24; m[32214] = 1.880606E-24; m[32224] = 1.880606E-24;
- m[33122] = 1.880606E-23; m[42112] = 6.582122E-24; m[42212] = 6.582122E-24;
- m[43122] = 2.194041E-24; m[53122] = 4.388081E-24; m[100111] = 1.645531E-24;
- m[100113] = 1.64553E-24; m[100211] = 1.645531E-24; m[100213] = 1.64553E-24;
- m[100221] = 1.196749E-23; m[100223] = 3.061452E-24; m[100313] = 2.837122E-24;
- m[100323] = 2.837122E-24; m[100331] = 4.459432E-25; m[100333] = 4.388081E-24;
- m[100441] = 4.701516E-23; m[100443] = 2.076379E-21; m[100553] = 2.056913E-20;
- m[200553] = 3.242425E-20; m[300553] = 3.210791E-23; m[9000111] = 8.776163E-24;
- m[9000211] = 8.776163E-24; m[9000443] = 8.227652E-24; m[9000553] = 5.983747E-24;
- m[9010111] = 3.164482E-24; m[9010211] = 3.164482E-24; m[9010221] = 9.403031E-24;
- m[9010443] = 8.438618E-24; m[9010553] = 8.3318E-24;
- m[9020221] = 8.093281E-23; m[9020443] = 1.061633E-23;
- m[9030221] = 6.038644E-24; m[9042413] = 2.07634E-21; m[9050225] = 1.394517E-24;
- m[9060225] = 3.291061E-24; m[9080225] = 4.388081E-24; m[9090225] = 2.056913E-24;
- m[9910445] = 2.07634E-21; m[9920443] = 2.07634E-21;
+ m[6] = 4.707703E-25;
+ m[11] = 1.E+16;
+ m[12] = 1.E+16;
+ m[13] = 2.197019E-06;
+ m[14] = 1.E+16;
+ m[15] = 2.906E-13;
+ m[16] = 1.E+16;
+ m[22] = 1.E+16;
+ m[23] = 2.637914E-25;
+ m[24] = 3.075758E-25;
+ m[25] = 9.4E-26;
+ m[35] = 9.4E-26;
+ m[36] = 9.4E-26;
+ m[37] = 9.4E-26;
+ m[84] = 3.335641E-13;
+ m[85] = 1.290893E-12;
+ m[111] = 8.4E-17;
+ m[113] = 4.405704E-24;
+ m[115] = 6.151516E-24;
+ m[117] = 4.088275E-24;
+ m[119] = 2.102914E-24;
+ m[130] = 5.116E-08;
+ m[150] = 1.525E-12;
+ m[211] = 2.6033E-08;
+ m[213] = 4.405704E-24;
+ m[215] = 6.151516E-24;
+ m[217] = 4.088275E-24;
+ m[219] = 2.102914E-24;
+ m[221] = 5.063171E-19;
+ m[223] = 7.752794E-23;
+ m[225] = 3.555982E-24;
+ m[227] = 3.91793E-24;
+ m[229] = 2.777267E-24;
+ m[310] = 8.953E-11;
+ m[313] = 1.308573E-23;
+ m[315] = 6.038644E-24;
+ m[317] = 4.139699E-24;
+ m[319] = 3.324304E-24;
+ m[321] = 1.238E-08;
+ m[323] = 1.295693E-23;
+ m[325] = 6.682357E-24;
+ m[327] = 4.139699E-24;
+ m[329] = 3.324304E-24;
+ m[331] = 3.210791E-21;
+ m[333] = 1.545099E-22;
+ m[335] = 9.016605E-24;
+ m[337] = 7.565657E-24;
+ m[350] = 1.407125E-12;
+ m[411] = 1.04E-12;
+ m[413] = 6.856377E-21;
+ m[415] = 1.778952E-23;
+ m[421] = 4.101E-13;
+ m[423] = 1.000003E-19;
+ m[425] = 1.530726E-23;
+ m[431] = 5.E-13;
+ m[433] = 1.000003E-19;
+ m[435] = 3.291061E-23;
+ m[441] = 2.465214E-23;
+ m[443] = 7.062363E-21;
+ m[445] = 3.242425E-22;
+ m[510] = 1.525E-12;
+ m[511] = 1.525E-12;
+ m[513] = 1.000019E-19;
+ m[515] = 1.31E-23;
+ m[521] = 1.638E-12;
+ m[523] = 1.000019E-19;
+ m[525] = 1.31E-23;
+ m[530] = 1.536875E-12;
+ m[531] = 1.472E-12;
+ m[533] = 1.E-19;
+ m[535] = 1.31E-23;
+ m[541] = 4.5E-13;
+ m[553] = 1.218911E-20;
+ m[1112] = 4.539394E-24;
+ m[1114] = 5.578069E-24;
+ m[1116] = 1.994582E-24;
+ m[1118] = 2.269697E-24;
+ m[1212] = 4.539394E-24;
+ m[1214] = 5.723584E-24;
+ m[1216] = 1.994582E-24;
+ m[1218] = 1.316424E-24;
+ m[2112] = 8.857E+02;
+ m[2114] = 5.578069E-24;
+ m[2116] = 4.388081E-24;
+ m[2118] = 2.269697E-24;
+ m[2122] = 4.539394E-24;
+ m[2124] = 5.723584E-24;
+ m[2126] = 1.994582E-24;
+ m[2128] = 1.316424E-24;
+ m[2212] = 1.E+16;
+ m[2214] = 5.578069E-24;
+ m[2216] = 4.388081E-24;
+ m[2218] = 2.269697E-24;
+ m[2222] = 4.539394E-24;
+ m[2224] = 5.578069E-24;
+ m[2226] = 1.994582E-24;
+ m[2228] = 2.269697E-24;
+ m[3112] = 1.479E-10;
+ m[3114] = 1.670589E-23;
+ m[3116] = 5.485102E-24;
+ m[3118] = 3.656734E-24;
+ m[3122] = 2.631E-10;
+ m[3124] = 4.219309E-23;
+ m[3126] = 8.227653E-24;
+ m[3128] = 3.291061E-24;
+ m[3212] = 7.4E-20;
+ m[3214] = 1.828367E-23;
+ m[3216] = 5.485102E-24;
+ m[3218] = 3.656734E-24;
+ m[3222] = 8.018E-11;
+ m[3224] = 1.838582E-23;
+ m[3226] = 5.485102E-24;
+ m[3228] = 3.656734E-24;
+ m[3312] = 1.639E-10;
+ m[3314] = 6.648608E-23;
+ m[3322] = 2.9E-10;
+ m[3324] = 7.233101E-23;
+ m[3334] = 8.21E-11;
+ m[4112] = 2.991874E-22;
+ m[4114] = 4.088274E-23;
+ m[4122] = 2.E-13;
+ m[4132] = 1.12E-13;
+ m[4212] = 3.999999E-22;
+ m[4214] = 3.291061E-22;
+ m[4222] = 2.951624E-22;
+ m[4224] = 4.417531E-23;
+ m[4232] = 4.42E-13;
+ m[4332] = 6.9E-14;
+ m[4412] = 3.335641E-13;
+ m[4422] = 3.335641E-13;
+ m[4432] = 3.335641E-13;
+ m[5112] = 1.E-19;
+ m[5122] = 1.38E-12;
+ m[5132] = 1.42E-12;
+ m[5142] = 1.290893E-12;
+ m[5212] = 1.E-19;
+ m[5222] = 1.E-19;
+ m[5232] = 1.42E-12;
+ m[5242] = 1.290893E-12;
+ m[5312] = 1.E-19;
+ m[5322] = 1.E-19;
+ m[5332] = 1.55E-12;
+ m[5342] = 1.290893E-12;
+ m[5442] = 1.290893E-12;
+ m[5512] = 1.290893E-12;
+ m[5522] = 1.290893E-12;
+ m[5532] = 1.290893E-12;
+ m[5542] = 1.290893E-12;
+ m[10111] = 2.48382E-24;
+ m[10113] = 4.635297E-24;
+ m[10115] = 2.54136E-24;
+ m[10211] = 2.48382E-24;
+ m[10213] = 4.635297E-24;
+ m[10215] = 2.54136E-24;
+ m[10223] = 1.828367E-24;
+ m[10225] = 3.636531E-24;
+ m[10311] = 2.437823E-24;
+ m[10313] = 7.313469E-24;
+ m[10315] = 3.538775E-24;
+ m[10321] = 2.437823E-24;
+ m[10323] = 7.313469E-24;
+ m[10325] = 3.538775E-24;
+ m[10331] = 4.804469E-24;
+ m[10411] = 4.38E-24;
+ m[10413] = 3.29E-23;
+ m[10421] = 4.38E-24;
+ m[10423] = 3.22653E-23;
+ m[10431] = 6.5821E-22;
+ m[10433] = 6.5821E-22;
+ m[10441] = 6.453061E-23;
+ m[10511] = 4.39E-24;
+ m[10513] = 1.65E-23;
+ m[10521] = 4.39E-24;
+ m[10523] = 1.65E-23;
+ m[10531] = 4.39E-24;
+ m[10533] = 1.65E-23;
+ m[11114] = 2.194041E-24;
+ m[11116] = 1.828367E-24;
+ m[11212] = 1.880606E-24;
+ m[11216] = 1.828367E-24;
+ m[12112] = 2.194041E-24;
+ m[12114] = 2.194041E-24;
+ m[12116] = 5.063171E-24;
+ m[12126] = 1.828367E-24;
+ m[12212] = 2.194041E-24;
+ m[12214] = 2.194041E-24;
+ m[12216] = 5.063171E-24;
+ m[12224] = 2.194041E-24;
+ m[12226] = 1.828367E-24;
+ m[13112] = 6.582122E-24;
+ m[13114] = 1.09702E-23;
+ m[13116] = 5.485102E-24;
+ m[13122] = 1.316424E-23;
+ m[13124] = 1.09702E-23;
+ m[13126] = 6.928549E-24;
+ m[13212] = 6.582122E-24;
+ m[13214] = 1.09702E-23;
+ m[13216] = 5.485102E-24;
+ m[13222] = 6.582122E-24;
+ m[13224] = 1.09702E-23;
+ m[13226] = 5.485102E-24;
+ m[13312] = 4.135667E-22;
+ m[13314] = 2.742551E-23;
+ m[13324] = 2.742551E-23;
+ m[14122] = 1.828367E-22;
+ m[20022] = 1.E+16;
+ m[20113] = 1.567172E-24;
+ m[20213] = 1.567172E-24;
+ m[20223] = 2.708692E-23;
+ m[20313] = 3.782829E-24;
+ m[20315] = 2.384827E-24;
+ m[20323] = 3.782829E-24;
+ m[20325] = 2.384827E-24;
+ m[20333] = 1.198929E-23;
+ m[20413] = 2.63E-24;
+ m[20423] = 2.63E-24;
+ m[20433] = 6.5821E-22;
+ m[20443] = 7.395643E-22;
+ m[20513] = 2.63E-24;
+ m[20523] = 2.63E-24;
+ m[20533] = 2.63E-24;
+ m[21112] = 2.632849E-24;
+ m[21114] = 3.291061E-24;
+ m[21212] = 2.632849E-24;
+ m[21214] = 6.582122E-24;
+ m[22112] = 4.388081E-24;
+ m[22114] = 3.291061E-24;
+ m[22122] = 2.632849E-24;
+ m[22124] = 6.582122E-24;
+ m[22212] = 4.388081E-24;
+ m[22214] = 3.291061E-24;
+ m[22222] = 2.632849E-24;
+ m[22224] = 3.291061E-24;
+ m[23112] = 7.313469E-24;
+ m[23114] = 2.991874E-24;
+ m[23122] = 4.388081E-24;
+ m[23124] = 6.582122E-24;
+ m[23126] = 3.291061E-24;
+ m[23212] = 7.313469E-24;
+ m[23214] = 2.991874E-24;
+ m[23222] = 7.313469E-24;
+ m[23224] = 2.991874E-24;
+ m[30113] = 2.632849E-24;
+ m[30213] = 2.632849E-24;
+ m[30221] = 1.880606E-24;
+ m[30223] = 2.089563E-24;
+ m[30313] = 2.056913E-24;
+ m[30323] = 2.056913E-24;
+ m[30443] = 2.419898E-23;
+ m[31114] = 1.880606E-24;
+ m[31214] = 3.291061E-24;
+ m[32112] = 3.989164E-24;
+ m[32114] = 1.880606E-24;
+ m[32124] = 3.291061E-24;
+ m[32212] = 3.989164E-24;
+ m[32214] = 1.880606E-24;
+ m[32224] = 1.880606E-24;
+ m[33122] = 1.880606E-23;
+ m[42112] = 6.582122E-24;
+ m[42212] = 6.582122E-24;
+ m[43122] = 2.194041E-24;
+ m[53122] = 4.388081E-24;
+ m[100111] = 1.645531E-24;
+ m[100113] = 1.64553E-24;
+ m[100211] = 1.645531E-24;
+ m[100213] = 1.64553E-24;
+ m[100221] = 1.196749E-23;
+ m[100223] = 3.061452E-24;
+ m[100313] = 2.837122E-24;
+ m[100323] = 2.837122E-24;
+ m[100331] = 4.459432E-25;
+ m[100333] = 4.388081E-24;
+ m[100441] = 4.701516E-23;
+ m[100443] = 2.076379E-21;
+ m[100553] = 2.056913E-20;
+ m[200553] = 3.242425E-20;
+ m[300553] = 3.210791E-23;
+ m[9000111] = 8.776163E-24;
+ m[9000211] = 8.776163E-24;
+ m[9000443] = 8.227652E-24;
+ m[9000553] = 5.983747E-24;
+ m[9010111] = 3.164482E-24;
+ m[9010211] = 3.164482E-24;
+ m[9010221] = 9.403031E-24;
+ m[9010443] = 8.438618E-24;
+ m[9010553] = 8.3318E-24;
+ m[9020221] = 8.093281E-23;
+ m[9020443] = 1.061633E-23;
+ m[9030221] = 6.038644E-24;
+ m[9042413] = 2.07634E-21;
+ m[9050225] = 1.394517E-24;
+ m[9060225] = 3.291061E-24;
+ m[9080225] = 4.388081E-24;
+ m[9090225] = 2.056913E-24;
+ m[9910445] = 2.07634E-21;
+ m[9920443] = 2.07634E-21;
return true;
}
/// @name Histograms
//@{
Histo1DPtr _h_K0s_pt_y_30; // histogram for 2.5 < y < 3.0 (d2sigma)
Histo1DPtr _h_K0s_pt_y_35; // histogram for 3.0 < y < 3.5 (d2sigma)
Histo1DPtr _h_K0s_pt_y_40; // histogram for 3.5 < y < 4.0 (d2sigma)
Histo1DPtr _h_K0s_pt_30; // histogram for 2.5 < y < 3.0 (sigma)
Histo1DPtr _h_K0s_pt_35; // histogram for 3.0 < y < 3.5 (sigma)
Histo1DPtr _h_K0s_pt_40; // histogram for 3.5 < y < 4.0 (sigma)
Histo1DPtr _h_K0s_pt_y_all; // histogram for 2.5 < y < 4.0 (d2sigma)
- double sumKs0_30; // Sum of weights 2.5 < y < 3.0
- double sumKs0_35; // Sum of weights 3.0 < y < 3.5
- double sumKs0_40; // Sum of weights 3.5 < y < 4.0
+ CounterPtr sumKs0_30; // Sum of weights 2.5 < y < 3.0
+ CounterPtr sumKs0_35; // Sum of weights 3.0 < y < 3.5
+ CounterPtr sumKs0_40; // Sum of weights 3.5 < y < 4.0
// Various counters mainly for debugging and comparisons between different generators
size_t sumKs0_badnull; // Nb of particles for which mother could not be identified
size_t sumKs0_badlft; // Nb of mesons with long lived mothers
size_t sumKs0_all; // Nb of all Ks0 generated
size_t sumKs0_outup; // Nb of mesons with y > 4.0
size_t sumKs0_outdwn; // Nb of mesons with y < 2.5
size_t sum_low_pt_loss; // Nb of mesons with very low pT (indicates when units are mixed-up)
size_t sum_high_pt_loss; // Nb of mesons with pT > 1.6 GeV/c
// Map between PDG id and particle lifetimes in seconds
std::map<int, double> partLftMap;
// Set of PDG Ids for stable particles (PDG Id <= 100 are considered stable)
- static const int stablePDGIds[205];
+ static const array<int,171> stablePDGIds;
//@}
};
// Actual initialization according to ISO C++ requirements
- const int LHCB_2010_S8758301::stablePDGIds[205] = {
+ const array<int,171> LHCB_2010_S8758301::stablePDGIds{{
311, 543, 545, 551, 555, 557, 1103, 2101, 2103, 2203, 3101, 3103, 3201, 3203, 3303,
4101, 4103, 4124, 4201, 4203, 4301, 4303, 4312, 4314, 4322, 4324, 4334, 4403, 4414,
4424, 4434, 4444, 5101, 5103, 5114, 5201, 5203, 5214, 5224, 5301, 5303, 5314, 5324,
5334, 5401, 5403, 5412, 5414, 5422, 5424, 5432, 5434, 5444, 5503, 5514, 5524, 5534,
5544, 5554, 10022, 10333, 10335, 10443, 10541, 10543, 10551, 10553, 10555, 11112,
12118, 12122, 12218, 12222, 13316, 13326, 20543, 20553, 20555, 23314, 23324, 30343,
30353, 30363, 30553, 33314, 33324, 41214, 42124, 52114, 52214, 100311, 100315, 100321,
100325, 100411, 100413, 100421, 100423, 100551, 100555, 100557, 110551, 110553, 110555,
120553, 120555, 130553, 200551, 200555, 210551, 210553, 220553, 1000001, 1000002,
1000003, 1000004, 1000005, 1000006, 1000011, 1000012, 1000013, 1000014, 1000015,
1000016, 1000021, 1000022, 1000023, 1000024, 1000025, 1000035, 1000037, 1000039,
2000001, 2000002, 2000003, 2000004, 2000005, 2000006, 2000011, 2000012, 2000013,
2000014, 2000015, 2000016, 3000111, 3000113, 3000211, 3000213, 3000221, 3000223,
3000331, 3100021, 3100111, 3100113, 3200111, 3200113, 3300113, 3400113, 4000001,
4000002, 4000011, 4000012, 5000039, 9000221, 9900012, 9900014, 9900016, 9900023,
- 9900024, 9900041, 9900042};
+ 9900024, 9900041, 9900042}};
// Hook for the plugin system
DECLARE_RIVET_PLUGIN(LHCB_2010_S8758301);
}
diff --git a/analyses/pluginLHCb/LHCB_2011_I917009.cc b/analyses/pluginLHCb/LHCB_2011_I917009.cc
--- a/analyses/pluginLHCb/LHCB_2011_I917009.cc
+++ b/analyses/pluginLHCb/LHCB_2011_I917009.cc
@@ -1,334 +1,333 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/UnstableFinalState.hh"
namespace Rivet {
class LHCB_2011_I917009 : public Analysis {
public:
/// @name Constructors etc.
//@{
/// Constructor
LHCB_2011_I917009()
: Analysis("LHCB_2011_I917009"),
rap_beam(0.0), pt_min(0.0),
pt1_edge(0.65), pt2_edge(1.0),
pt3_edge(2.5), rap_min(2.),
rap_max(0.0), dsShift(0)
{ }
//@}
public:
/// @name Analysis methods
//@{
/// Book histograms and initialise projections before the run
void init() {
int y_nbins = 4;
fillMap(partLftMap);
if (fuzzyEquals(sqrtS(), 0.9*TeV)) {
rap_beam = 6.87;
rap_max = 4.;
pt_min = 0.25;
} else if (fuzzyEquals(sqrtS(), 7*TeV)) {
rap_beam = 8.92;
rap_max = 4.5;
pt_min = 0.15;
y_nbins = 5;
dsShift = 8;
} else {
MSG_ERROR("Incompatible beam energy!");
}
// Create the sets of temporary histograms that will be used to make the ratios in the finalize()
for (size_t i = 0; i < 12; ++i) book(_tmphistos[i], "TMP/"+to_str(i), y_nbins, rap_min, rap_max);
for (size_t i = 12; i < 15; ++i) book(_tmphistos[i], "TMP/"+to_str(i), refData(dsShift+5, 1, 1));
for (size_t i = 15; i < 18; ++i) book(_tmphistos[i], "TMP/"+to_str(i), y_nbins, rap_beam - rap_max, rap_beam - rap_min);
int dsId = dsShift + 1;
for (size_t j = 0; j < 3; ++j) {
book(s1, dsId, 1, j+1);
book(s2, dsId+1, 1, j+1);
}
dsId += 2;
for (size_t j = 3; j < 6; ++j) {
book(s3, dsId, 1, 1);
dsId += 1;
book(s4, dsId, 1, 1);
dsId += 1;
}
declare(UnstableFinalState(), "UFS");
}
/// Perform the per-event analysis
void analyze(const Event& event) {
- const double weight = 1.0;
const UnstableFinalState& ufs = apply<UnstableFinalState>(event, "UFS");
double ancestor_lftsum = 0.0;
double y, pT;
int id;
int partIdx = -1;
foreach (const Particle& p, ufs.particles()) {
id = p.pid();
// continue if particle not a K0s nor (anti-)Lambda
if ( (id == 310) || (id == -310) ) {
partIdx = 2;
} else if ( id == 3122 ) {
partIdx = 1;
} else if ( id == -3122 ) {
partIdx = 0;
} else {
continue;
}
ancestor_lftsum = getMotherLifeTimeSum(p);
// Lifetime cut: ctau sum of all particle ancestors < 10^-9 m according to the paper (see eq. 5)
const double MAX_CTAU = 1.0E-9; // [m]
if ( (ancestor_lftsum < 0.0) || (ancestor_lftsum > MAX_CTAU) ) continue;
const FourMomentum& qmom = p.momentum();
y = log((qmom.E() + qmom.pz())/(qmom.E() - qmom.pz()))/2.;
// skip this particle if it has too high or too low rapidity (extremely rare cases when E = +- pz)
if ( std::isnan(y) || std::isinf(y) ) continue;
y = fabs(y);
if (!inRange(y, rap_min, rap_max)) continue;
pT = sqrt((qmom.px() * qmom.px()) + (qmom.py() * qmom.py()));
if (!inRange(pT, pt_min, pt3_edge)) continue;
// Filling corresponding temporary histograms for pT intervals
- if (inRange(pT, pt_min, pt1_edge)) _tmphistos[partIdx*3]->fill(y, weight);
- if (inRange(pT, pt1_edge, pt2_edge)) _tmphistos[partIdx*3+1]->fill(y, weight);
- if (inRange(pT, pt2_edge, pt3_edge)) _tmphistos[partIdx*3+2]->fill(y, weight);
+ if (inRange(pT, pt_min, pt1_edge)) _tmphistos[partIdx*3]->fill(y);
+ if (inRange(pT, pt1_edge, pt2_edge)) _tmphistos[partIdx*3+1]->fill(y);
+ if (inRange(pT, pt2_edge, pt3_edge)) _tmphistos[partIdx*3+2]->fill(y);
// Fill histo in rapidity for whole pT interval
- _tmphistos[partIdx+9]->fill(y, weight);
+ _tmphistos[partIdx+9]->fill(y);
// Fill histo in pT for whole rapidity interval
- _tmphistos[partIdx+12]->fill(pT, weight);
+ _tmphistos[partIdx+12]->fill(pT);
// Fill histo in rapidity loss for whole pT interval
- _tmphistos[partIdx+15]->fill(rap_beam - y, weight);
+ _tmphistos[partIdx+15]->fill(rap_beam - y);
}
}
// Generate the ratio histograms
void finalize() {
int dsId = dsShift + 1;
for (size_t j = 0; j < 3; ++j) {
/// @todo Compactify to two one-liners
divide(_tmphistos[j], _tmphistos[3+j], s1);
divide(_tmphistos[j], _tmphistos[6+j], s2);
}
dsId += 2;
for (size_t j = 3; j < 6; ++j) {
/// @todo Compactify to two one-liners
divide(_tmphistos[3*j], _tmphistos[3*j+1], s3);
dsId += 1;
divide(_tmphistos[3*j], _tmphistos[3*j+2], s4);
dsId += 1;
}
}
//@}
private:
// Get particle lifetime from hardcoded data
double getLifeTime(int pid) {
double lft = -1.0;
if (pid < 0) pid = - pid;
// Correct Pythia6 PIDs for f0(980), f0(1370) mesons
if (pid == 10331) pid = 30221;
if (pid == 10221) pid = 9010221;
map<int, double>::iterator pPartLft = partLftMap.find(pid);
// search stable particle list
if (pPartLft == partLftMap.end()) {
if (pid <= 100) return 0.0;
for (size_t i=0; i < sizeof(stablePDGIds)/sizeof(unsigned int); i++) {
if (pid == stablePDGIds[i]) { lft = 0.0; break; }
}
} else {
lft = (*pPartLft).second;
}
if (lft < 0.0 && PID::isHadron(pid)) {
MSG_ERROR("Could not determine lifetime for particle with PID " << pid
<< "... This V^0 will be considered unprompt!");
}
return lft;
}
// Data members like post-cuts event weight counters go here
const double getMotherLifeTimeSum(const Particle& p) {
if (p.genParticle() == NULL) return -1.;
double lftSum = 0.;
double plft = 0.;
const GenParticle* part = p.genParticle();
const GenVertex* ivtx = part->production_vertex();
while (ivtx) {
if (ivtx->particles_in_size() < 1) { lftSum = -1.; break; };
const GenVertex::particles_in_const_iterator iPart_invtx = ivtx->particles_in_const_begin();
part = (*iPart_invtx);
if ( !(part) ) { lftSum = -1.; break; };
ivtx = part->production_vertex();
if ( (part->pdg_id() == 2212) || !(ivtx) ) break; //reached beam
plft = getLifeTime(part->pdg_id());
if (plft < 0.) { lftSum = -1.; break; };
lftSum += plft;
};
return (lftSum * c_light);
}
/// @name Private variables
//@{
// The rapidity of the beam according to the selected beam energy
double rap_beam;
// The edges of the intervals of transverse momentum
double pt_min, pt1_edge, pt2_edge, pt3_edge;
// The limits of the rapidity window
double rap_min;
double rap_max;
// Indicates which set of histograms will be output to yoda file (according to beam energy)
int dsShift;
// Map between PDG id and particle lifetimes in seconds
std::map<int, double> partLftMap;
// Set of PDG Ids for stable particles (PDG Id <= 100 are considered stable)
static const int stablePDGIds[205];
//@}
/// @name Helper histograms
//@{
/// Histograms are defined in the following order: anti-Lambda, Lambda and K0s.
/// First 3 suites of 3 histograms correspond to each particle in bins of y for the 3 pT intervals. (9 histos)
/// Next 3 histograms contain the particles in y bins for the whole pT interval (3 histos)
/// Next 3 histograms contain the particles in y_loss bins for the whole pT interval (3 histos)
/// Last 3 histograms contain the particles in pT bins for the whole rapidity (y) interval (3 histos)
Histo1DPtr _tmphistos[18];
Scatter2DPtr s1,s2,s3,s4;
//@}
// Fill the PDG Id to Lifetime[seconds] map
// Data was extracted from LHCb Particle Table through LHCb::ParticlePropertySvc
bool fillMap(map<int, double>& m) {
m[6] = 4.707703E-25; m[11] = 1.E+16; m[12] = 1.E+16;
m[13] = 2.197019E-06; m[14] = 1.E+16; m[15] = 2.906E-13; m[16] = 1.E+16; m[22] = 1.E+16;
m[23] = 2.637914E-25; m[24] = 3.075758E-25; m[25] = 9.4E-26; m[35] = 9.4E-26;
m[36] = 9.4E-26; m[37] = 9.4E-26; m[84] = 3.335641E-13; m[85] = 1.290893E-12;
m[111] = 8.4E-17; m[113] = 4.405704E-24; m[115] = 6.151516E-24; m[117] = 4.088275E-24;
m[119] = 2.102914E-24; m[130] = 5.116E-08; m[150] = 1.525E-12; m[211] = 2.6033E-08;
m[213] = 4.405704E-24; m[215] = 6.151516E-24; m[217] = 4.088275E-24; m[219] = 2.102914E-24;
m[221] = 5.063171E-19; m[223] = 7.752794E-23; m[225] = 3.555982E-24; m[227] = 3.91793E-24;
m[229] = 2.777267E-24; m[310] = 8.953E-11; m[313] = 1.308573E-23; m[315] = 6.038644E-24;
m[317] = 4.139699E-24; m[319] = 3.324304E-24; m[321] = 1.238E-08; m[323] = 1.295693E-23;
m[325] = 6.682357E-24; m[327] = 4.139699E-24; m[329] = 3.324304E-24; m[331] = 3.210791E-21;
m[333] = 1.545099E-22; m[335] = 9.016605E-24; m[337] = 7.565657E-24; m[350] = 1.407125E-12;
m[411] = 1.04E-12; m[413] = 6.856377E-21; m[415] = 1.778952E-23; m[421] = 4.101E-13;
m[423] = 1.000003E-19; m[425] = 1.530726E-23; m[431] = 5.E-13; m[433] = 1.000003E-19;
m[435] = 3.291061E-23; m[441] = 2.465214E-23; m[443] = 7.062363E-21; m[445] = 3.242425E-22;
m[510] = 1.525E-12; m[511] = 1.525E-12; m[513] = 1.000019E-19; m[515] = 1.31E-23;
m[521] = 1.638E-12; m[523] = 1.000019E-19; m[525] = 1.31E-23; m[530] = 1.536875E-12;
m[531] = 1.472E-12; m[533] = 1.E-19; m[535] = 1.31E-23; m[541] = 4.5E-13;
m[553] = 1.218911E-20; m[1112] = 4.539394E-24; m[1114] = 5.578069E-24; m[1116] = 1.994582E-24;
m[1118] = 2.269697E-24; m[1212] = 4.539394E-24; m[1214] = 5.723584E-24; m[1216] = 1.994582E-24;
m[1218] = 1.316424E-24; m[2112] = 8.857E+02; m[2114] = 5.578069E-24; m[2116] = 4.388081E-24;
m[2118] = 2.269697E-24; m[2122] = 4.539394E-24; m[2124] = 5.723584E-24; m[2126] = 1.994582E-24;
m[2128] = 1.316424E-24; m[2212] = 1.E+16; m[2214] = 5.578069E-24; m[2216] = 4.388081E-24;
m[2218] = 2.269697E-24; m[2222] = 4.539394E-24; m[2224] = 5.578069E-24; m[2226] = 1.994582E-24;
m[2228] = 2.269697E-24; m[3112] = 1.479E-10; m[3114] = 1.670589E-23; m[3116] = 5.485102E-24;
m[3118] = 3.656734E-24; m[3122] = 2.631E-10; m[3124] = 4.219309E-23; m[3126] = 8.227653E-24;
m[3128] = 3.291061E-24; m[3212] = 7.4E-20; m[3214] = 1.828367E-23; m[3216] = 5.485102E-24;
m[3218] = 3.656734E-24; m[3222] = 8.018E-11; m[3224] = 1.838582E-23; m[3226] = 5.485102E-24;
m[3228] = 3.656734E-24; m[3312] = 1.639E-10; m[3314] = 6.648608E-23; m[3322] = 2.9E-10;
m[3324] = 7.233101E-23; m[3334] = 8.21E-11; m[4112] = 2.991874E-22; m[4114] = 4.088274E-23;
m[4122] = 2.E-13; m[4132] = 1.12E-13; m[4212] = 3.999999E-22; m[4214] = 3.291061E-22;
m[4222] = 2.951624E-22; m[4224] = 4.417531E-23; m[4232] = 4.42E-13; m[4332] = 6.9E-14;
m[4412] = 3.335641E-13; m[4422] = 3.335641E-13; m[4432] = 3.335641E-13; m[5112] = 1.E-19;
m[5122] = 1.38E-12; m[5132] = 1.42E-12; m[5142] = 1.290893E-12; m[5212] = 1.E-19;
m[5222] = 1.E-19; m[5232] = 1.42E-12; m[5242] = 1.290893E-12; m[5312] = 1.E-19;
m[5322] = 1.E-19; m[5332] = 1.55E-12; m[5342] = 1.290893E-12; m[5442] = 1.290893E-12;
m[5512] = 1.290893E-12; m[5522] = 1.290893E-12; m[5532] = 1.290893E-12; m[5542] = 1.290893E-12;
m[10111] = 2.48382E-24; m[10113] = 4.635297E-24; m[10115] = 2.54136E-24; m[10211] = 2.48382E-24;
m[10213] = 4.635297E-24; m[10215] = 2.54136E-24; m[10223] = 1.828367E-24; m[10225] = 3.636531E-24;
m[10311] = 2.437823E-24; m[10313] = 7.313469E-24; m[10315] = 3.538775E-24;
m[10321] = 2.437823E-24; m[10323] = 7.313469E-24; m[10325] = 3.538775E-24;
m[10331] = 4.804469E-24; m[10411] = 4.38E-24; m[10413] = 3.29E-23; m[10421] = 4.38E-24;
m[10423] = 3.22653E-23; m[10431] = 6.5821E-22; m[10433] = 6.5821E-22; m[10441] = 6.453061E-23;
m[10511] = 4.39E-24; m[10513] = 1.65E-23; m[10521] = 4.39E-24; m[10523] = 1.65E-23;
m[10531] = 4.39E-24; m[10533] = 1.65E-23; m[11114] = 2.194041E-24; m[11116] = 1.828367E-24;
m[11212] = 1.880606E-24; m[11216] = 1.828367E-24; m[12112] = 2.194041E-24;
m[12114] = 2.194041E-24; m[12116] = 5.063171E-24; m[12126] = 1.828367E-24;
m[12212] = 2.194041E-24; m[12214] = 2.194041E-24; m[12216] = 5.063171E-24;
m[12224] = 2.194041E-24; m[12226] = 1.828367E-24; m[13112] = 6.582122E-24; m[13114] = 1.09702E-23;
m[13116] = 5.485102E-24; m[13122] = 1.316424E-23; m[13124] = 1.09702E-23; m[13126] = 6.928549E-24;
m[13212] = 6.582122E-24; m[13214] = 1.09702E-23; m[13216] = 5.485102E-24; m[13222] = 6.582122E-24;
m[13224] = 1.09702E-23; m[13226] = 5.485102E-24; m[13314] = 2.742551E-23;
m[13324] = 2.742551E-23; m[14122] = 1.828367E-22; m[20022] = 1.E+16; m[20113] = 1.567172E-24;
m[20213] = 1.567172E-24; m[20223] = 2.708692E-23; m[20313] = 3.782829E-24;
m[20315] = 2.384827E-24; m[20323] = 3.782829E-24; m[20325] = 2.384827E-24;
m[20333] = 1.198929E-23; m[20413] = 2.63E-24; m[20423] = 2.63E-24; m[20433] = 6.5821E-22;
m[20443] = 7.395643E-22; m[20513] = 2.63E-24; m[20523] = 2.63E-24; m[20533] = 2.63E-24;
m[21112] = 2.632849E-24; m[21114] = 3.291061E-24; m[21212] = 2.632849E-24;
m[21214] = 6.582122E-24; m[22112] = 4.388081E-24; m[22114] = 3.291061E-24;
m[22122] = 2.632849E-24; m[22124] = 6.582122E-24; m[22212] = 4.388081E-24;
m[22214] = 3.291061E-24; m[22222] = 2.632849E-24; m[22224] = 3.291061E-24;
m[23112] = 7.313469E-24; m[23114] = 2.991874E-24; m[23122] = 4.388081E-24;
m[23124] = 6.582122E-24; m[23126] = 3.291061E-24; m[23212] = 7.313469E-24;
m[23214] = 2.991874E-24; m[23222] = 7.313469E-24; m[23224] = 2.991874E-24;
m[30113] = 2.632849E-24; m[30213] = 2.632849E-24; m[30221] = 1.880606E-24;
m[30223] = 2.089563E-24; m[30313] = 2.056913E-24; m[30323] = 2.056913E-24;
m[30443] = 2.419898E-23; m[31114] = 1.880606E-24; m[31214] = 3.291061E-24;
m[32112] = 3.989164E-24; m[32114] = 1.880606E-24; m[32124] = 3.291061E-24;
m[32212] = 3.989164E-24; m[32214] = 1.880606E-24; m[32224] = 1.880606E-24;
m[33122] = 1.880606E-23; m[42112] = 6.582122E-24; m[42212] = 6.582122E-24;
m[43122] = 2.194041E-24; m[53122] = 4.388081E-24; m[100111] = 1.645531E-24;
m[100113] = 1.64553E-24; m[100211] = 1.645531E-24; m[100213] = 1.64553E-24;
m[100221] = 1.196749E-23; m[100223] = 3.061452E-24; m[100313] = 2.837122E-24;
m[100323] = 2.837122E-24; m[100331] = 4.459432E-25; m[100333] = 4.388081E-24;
m[100441] = 4.701516E-23; m[100443] = 2.076379E-21; m[100553] = 2.056913E-20;
m[200553] = 3.242425E-20; m[300553] = 3.210791E-23; m[9000111] = 8.776163E-24;
m[9000211] = 8.776163E-24; m[9000443] = 8.227652E-24; m[9000553] = 5.983747E-24;
m[9010111] = 3.164482E-24; m[9010211] = 3.164482E-24; m[9010221] = 9.403031E-24;
m[9010443] = 8.438618E-24; m[9010553] = 8.3318E-24; m[9020443] = 1.061633E-23;
m[9030221] = 6.038644E-24; m[9042413] = 2.07634E-21; m[9050225] = 1.394517E-24;
m[9060225] = 3.291061E-24; m[9080225] = 4.388081E-24; m[9090225] = 2.056913E-24;
m[9910445] = 2.07634E-21; m[9920443] = 2.07634E-21;
return true;
}
};
const int LHCB_2011_I917009::stablePDGIds[205] = {
311, 543, 545, 551, 555, 557, 1103, 2101, 2103, 2203, 3101, 3103, 3201, 3203, 3303,
4101, 4103, 4124, 4201, 4203, 4301, 4303, 4312, 4314, 4322, 4324, 4334, 4403, 4414,
4424, 4434, 4444, 5101, 5103, 5114, 5201, 5203, 5214, 5224, 5301, 5303, 5314, 5324,
5334, 5401, 5403, 5412, 5414, 5422, 5424, 5432, 5434, 5444, 5503, 5514, 5524, 5534,
5544, 5554, 10022, 10333, 10335, 10443, 10541, 10543, 10551, 10553, 10555, 11112,
12118, 12122, 12218, 12222, 13316, 13326, 20543, 20553, 20555, 23314, 23324, 30343,
30353, 30363, 30553, 33314, 33324, 41214, 42124, 52114, 52214, 100311, 100315, 100321,
100325, 100411, 100413, 100421, 100423, 100551, 100555, 100557, 110551, 110553, 110555,
120553, 120555, 130553, 200551, 200555, 210551, 210553, 220553, 1000001, 1000002,
1000003, 1000004, 1000005, 1000006, 1000011, 1000012, 1000013, 1000014, 1000015,
1000016, 1000021, 1000022, 1000023, 1000024, 1000025, 1000035, 1000037, 1000039,
2000001, 2000002, 2000003, 2000004, 2000005, 2000006, 2000011, 2000012, 2000013,
2000014, 2000015, 2000016, 3000111, 3000113, 3000211, 3000213, 3000221, 3000223,
3000331, 3100021, 3100111, 3100113, 3200111, 3200113, 3300113, 3400113, 4000001,
4000002, 4000011, 4000012, 5000039, 9000221, 9900012, 9900014, 9900016, 9900023,
9900024, 9900041, 9900042 };
// Hook for the plugin system
DECLARE_RIVET_PLUGIN(LHCB_2011_I917009);
}
diff --git a/analyses/pluginLHCb/LHCB_2012_I1119400.cc b/analyses/pluginLHCb/LHCB_2012_I1119400.cc
--- a/analyses/pluginLHCb/LHCB_2012_I1119400.cc
+++ b/analyses/pluginLHCb/LHCB_2012_I1119400.cc
@@ -1,357 +1,356 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/ChargedFinalState.hh"
namespace Rivet {
class LHCB_2012_I1119400 : public Analysis {
public:
/// @name Constructors etc.
//@{
/// Constructor
LHCB_2012_I1119400() : Analysis("LHCB_2012_I1119400"),
_p_min(5.0),
_pt_min(0.0),_pt1_edge(0.8), _pt2_edge(1.2),
//_eta_nbins(4),
_eta_min(2.5),
_eta_max(4.5)
{ }
//@}
public:
/// @name Analysis methods
//@{
/// Book histograms and initialise projections before the run
void init() {
fillMap(_partLftMap);
int id_shift = 0;
if (fuzzyEquals(sqrtS(), 7*TeV)) id_shift = 1;
// define ratios if second pdgid in pair is -1, it means that is a antiparticle/particle ratio
_ratiotype["pbarp"] = make_pair(2212, -1);
_ratiotype["kminuskplus"] = make_pair(321, -1);
_ratiotype["piminuspiplus"] = make_pair(211, -1);
_ratiotype["ppi"] = make_pair(2212, 211);
_ratiotype["kpi"] = make_pair(321, 211);
_ratiotype["pk"] = make_pair(2212, 321);
std::map<string, int > _hepdataid;
_hepdataid["pbarp"] = 1 + id_shift;
_hepdataid["kminuskplus"] = 3 + id_shift;
_hepdataid["piminuspiplus"] = 5 + id_shift;
_hepdataid["ppi"] = 7 + id_shift;
_hepdataid["kpi"] = 9 + id_shift;
_hepdataid["pk"] = 11 + id_shift;
std::map<std::string, std::pair<int, int> >::iterator it;
// booking histograms
for (it=_ratiotype.begin(); it!=_ratiotype.end(); it++) {
book(_h_ratio_lowpt [it->first], _hepdataid[it->first], 1, 1);
book(_h_ratio_midpt [it->first], _hepdataid[it->first], 1, 2);
book(_h_ratio_highpt[it->first], _hepdataid[it->first], 1, 3);
book(_h_num_lowpt [it->first], "TMP/num_l_"+it->first,refData(_hepdataid[it->first], 1, 1));
book(_h_num_midpt [it->first], "TMP/num_m_"+it->first,refData(_hepdataid[it->first], 1, 2));
book(_h_num_highpt [it->first], "TMP/num_h_"+it->first,refData(_hepdataid[it->first], 1, 3));
book(_h_den_lowpt [it->first], "TMP/den_l_"+it->first,refData(_hepdataid[it->first], 1, 1));
book(_h_den_midpt [it->first], "TMP/den_m_"+it->first,refData(_hepdataid[it->first], 1, 2));
book(_h_den_highpt [it->first], "TMP/den_h_"+it->first,refData(_hepdataid[it->first], 1, 3));
}
declare(ChargedFinalState(_eta_min, _eta_max, _pt_min*GeV), "CFS");
}
// Perform the per-event analysis
void analyze(const Event& event) {
- const double weight = 1.0;
const ChargedFinalState& cfs = apply<ChargedFinalState>(event, "CFS");
foreach (const Particle& p, cfs.particles()) {
int id = p.pid();
// continue if particle not a proton, a kaon or a pion
if ( !( (abs(id) == 211) || (abs(id) == 321) || (abs(id) == 2212))) {
continue;
}
// cut in momentum
const FourMomentum& qmom = p.momentum();
if (qmom.p3().mod() < _p_min) continue;
// Lifetime cut: ctau sum of all particle ancestors < 10^-9 m according to the paper (see eq. 5)
const double MAX_CTAU = 1.0e-9; // [m]
double ancestor_lftsum = getMotherLifeTimeSum(p);
if ( (ancestor_lftsum < 0.0) || (ancestor_lftsum > MAX_CTAU) ) continue;
double eta = qmom.eta();
double pT = qmom.pT();
std::map<std::string, std::pair<int, int> >::iterator it;
for (it=_ratiotype.begin(); it!=_ratiotype.end(); it++) {
// check what type of ratio is
if ((it->second.second)==-1) {
// check ptbin
if (pT < _pt1_edge) {
// filling histos for numerator and denominator
- if (id == -abs(it->second.first)) _h_num_lowpt[it->first]->fill(eta, weight);
- if (id == abs(it->second.first)) _h_den_lowpt[it->first]->fill(eta, weight);
+ if (id == -abs(it->second.first)) _h_num_lowpt[it->first]->fill(eta);
+ if (id == abs(it->second.first)) _h_den_lowpt[it->first]->fill(eta);
}
else if (pT < _pt2_edge) {
// filling histos for numerator and denominator
- if (id == -abs(it->second.first)) _h_num_midpt[it->first]->fill(eta, weight);
- if (id == abs(it->second.first)) _h_den_midpt[it->first]->fill(eta, weight);
+ if (id == -abs(it->second.first)) _h_num_midpt[it->first]->fill(eta);
+ if (id == abs(it->second.first)) _h_den_midpt[it->first]->fill(eta);
}
else {
// filling histos for numerator and denominator
- if (id == -abs(it->second.first)) _h_num_highpt[it->first]->fill(eta, weight);
- if (id == abs(it->second.first)) _h_den_highpt[it->first]->fill(eta, weight);
+ if (id == -abs(it->second.first)) _h_num_highpt[it->first]->fill(eta);
+ if (id == abs(it->second.first)) _h_den_highpt[it->first]->fill(eta);
}
}
else {
// check what type of ratio is
if (pT < _pt1_edge) {
// filling histos for numerator and denominator
- if (abs(id) == abs(it->second.first)) _h_num_lowpt[it->first]->fill(eta, weight);
- if (abs(id) == abs(it->second.second)) _h_den_lowpt[it->first]->fill(eta, weight);
+ if (abs(id) == abs(it->second.first)) _h_num_lowpt[it->first]->fill(eta);
+ if (abs(id) == abs(it->second.second)) _h_den_lowpt[it->first]->fill(eta);
}
else if (pT < _pt2_edge) {
// filling histos for numerator and denominator
- if (abs(id) == abs(it->second.first)) _h_num_midpt[it->first]->fill(eta, weight);
- if (abs(id) == abs(it->second.second)) _h_den_midpt[it->first]->fill(eta, weight);
+ if (abs(id) == abs(it->second.first)) _h_num_midpt[it->first]->fill(eta);
+ if (abs(id) == abs(it->second.second)) _h_den_midpt[it->first]->fill(eta);
}
else {
// filling histos for numerator and denominator
- if (abs(id) == abs(it->second.first)) _h_num_highpt[it->first]->fill(eta, weight);
- if (abs(id) == abs(it->second.second)) _h_den_highpt[it->first]->fill(eta, weight);
+ if (abs(id) == abs(it->second.first)) _h_num_highpt[it->first]->fill(eta);
+ if (abs(id) == abs(it->second.second)) _h_den_highpt[it->first]->fill(eta);
}
}
}
}
}
// Generate the ratio histograms
void finalize() {
std::map<std::string, std::pair<int, int> >::iterator it;
// booking histograms
for (it=_ratiotype.begin(); it!=_ratiotype.end(); it++) {
divide(_h_num_lowpt[it->first], _h_den_lowpt[it->first], _h_ratio_lowpt[it->first]);
divide(_h_num_midpt[it->first], _h_den_midpt[it->first], _h_ratio_midpt[it->first]);
divide(_h_num_highpt[it->first], _h_den_highpt[it->first], _h_ratio_highpt[it->first]);
}
}
//@}
private:
// Get particle lifetime from hardcoded data
double getLifeTime(int pid) {
pid = abs(pid);
double lft = -1.0;
map<int, double>::iterator pPartLft = _partLftMap.find(pid);
// search stable particle list
if (pPartLft == _partLftMap.end()) {
if (pid <= 100) return 0.0;
for (size_t i=0; i < sizeof(_stablePDGIds)/sizeof(unsigned int); i++) {
if (pid == _stablePDGIds[i]) {
lft = 0.0;
break;
}
}
} else {
lft = (*pPartLft).second;
}
if (lft < 0.0 && PID::isHadron(pid)) {
MSG_WARNING("Lifetime map imcomplete --- " << pid
<< "... assume zero lifetime");
lft = 0.0;
}
return lft;
}
// Data members like post-cuts event weight counters go here
const double getMotherLifeTimeSum(const Particle& p) {
if (p.genParticle() == NULL) return -1.;
double lftSum = 0.;
double plft = 0.;
const GenParticle* part = p.genParticle();
const GenVertex* ivtx = part->production_vertex();
while(ivtx)
{
if (ivtx->particles_in_size() < 1) { lftSum = -1.; break; };
const GenVertex::particles_in_const_iterator iPart_invtx = ivtx->particles_in_const_begin();
part = (*iPart_invtx);
if ( !(part) ) { lftSum = -1.; break; };
ivtx = part->production_vertex();
if ( (part->pdg_id() == 2212) || !(ivtx) ) break; // reached beam
plft = getLifeTime(part->pdg_id());
if (plft < 0.) { lftSum = -1.; break; };
lftSum += plft;
};
return (lftSum * c_light);
}
/// @name Private variables
// Momentum threshold
double _p_min;
// The edges of the intervals of transversal momentum
double _pt_min;
double _pt1_edge;
double _pt2_edge;
// The limits of the pseudorapidity window
//int _eta_nbins;
double _eta_min;
double _eta_max;
// Map between PDG id and particle lifetimes in seconds
std::map<int, double> _partLftMap;
// Set of PDG Ids for stable particles (PDG Id <= 100 are considered stable)
static const int _stablePDGIds[205];
// Define histograms
// ratio
std::map<std::string, Scatter2DPtr > _h_ratio_lowpt;
std::map<std::string, Scatter2DPtr > _h_ratio_midpt;
std::map<std::string, Scatter2DPtr > _h_ratio_highpt;
// numerator
std::map<std::string, Histo1DPtr > _h_num_lowpt;
std::map<std::string, Histo1DPtr > _h_num_midpt;
std::map<std::string, Histo1DPtr > _h_num_highpt;
// denominator
std::map<std::string, Histo1DPtr > _h_den_lowpt;
std::map<std::string, Histo1DPtr > _h_den_midpt;
std::map<std::string, Histo1DPtr > _h_den_highpt;
// Map of ratios and IDs of numerator and denominator
std::map<string, pair<int,int> > _ratiotype;
// Fill the PDG Id to Lifetime[seconds] map
// Data was extracted from LHCb Particle Table through LHCb::ParticlePropertySvc
bool fillMap(map<int, double> &m) {
m[6] = 4.707703E-25; m[11] = 1.E+16; m[12] = 1.E+16;
m[13] = 2.197019E-06; m[14] = 1.E+16; m[15] = 2.906E-13; m[16] = 1.E+16; m[22] = 1.E+16;
m[23] = 2.637914E-25; m[24] = 3.075758E-25; m[25] = 9.4E-26; m[35] = 9.4E-26;
m[36] = 9.4E-26; m[37] = 9.4E-26; m[84] = 3.335641E-13; m[85] = 1.290893E-12;
m[111] = 8.4E-17; m[113] = 4.405704E-24; m[115] = 6.151516E-24; m[117] = 4.088275E-24;
m[119] = 2.102914E-24; m[130] = 5.116E-08; m[150] = 1.525E-12; m[211] = 2.6033E-08;
m[213] = 4.405704E-24; m[215] = 6.151516E-24; m[217] = 4.088275E-24; m[219] = 2.102914E-24;
m[221] = 5.063171E-19; m[223] = 7.752794E-23; m[225] = 3.555982E-24; m[227] = 3.91793E-24;
m[229] = 2.777267E-24; m[310] = 8.953E-11; m[313] = 1.308573E-23; m[315] = 6.038644E-24;
m[317] = 4.139699E-24; m[319] = 3.324304E-24; m[321] = 1.238E-08; m[323] = 1.295693E-23;
m[325] = 6.682357E-24; m[327] = 4.139699E-24; m[329] = 3.324304E-24; m[331] = 3.210791E-21;
m[333] = 1.545099E-22; m[335] = 9.016605E-24; m[337] = 7.565657E-24; m[350] = 1.407125E-12;
m[411] = 1.04E-12; m[413] = 6.856377E-21; m[415] = 1.778952E-23; m[421] = 4.101E-13;
m[423] = 1.000003E-19; m[425] = 1.530726E-23; m[431] = 5.E-13; m[433] = 1.000003E-19;
m[435] = 3.291061E-23; m[441] = 2.465214E-23; m[443] = 7.062363E-21; m[445] = 3.242425E-22;
m[510] = 1.525E-12; m[511] = 1.525E-12; m[513] = 1.000019E-19; m[515] = 1.31E-23;
m[521] = 1.638E-12; m[523] = 1.000019E-19; m[525] = 1.31E-23; m[530] = 1.536875E-12;
m[531] = 1.472E-12; m[533] = 1.E-19; m[535] = 1.31E-23; m[541] = 4.5E-13;
m[553] = 1.218911E-20; m[1112] = 4.539394E-24; m[1114] = 5.578069E-24; m[1116] = 1.994582E-24;
m[1118] = 2.269697E-24; m[1212] = 4.539394E-24; m[1214] = 5.723584E-24; m[1216] = 1.994582E-24;
m[1218] = 1.316424E-24; m[2112] = 8.857E+02; m[2114] = 5.578069E-24; m[2116] = 4.388081E-24;
m[2118] = 2.269697E-24; m[2122] = 4.539394E-24; m[2124] = 5.723584E-24; m[2126] = 1.994582E-24;
m[2128] = 1.316424E-24; m[2212] = 1.E+16; m[2214] = 5.578069E-24; m[2216] = 4.388081E-24;
m[2218] = 2.269697E-24; m[2222] = 4.539394E-24; m[2224] = 5.578069E-24; m[2226] = 1.994582E-24;
m[2228] = 2.269697E-24; m[3112] = 1.479E-10; m[3114] = 1.670589E-23; m[3116] = 5.485102E-24;
m[3118] = 3.656734E-24; m[3122] = 2.631E-10; m[3124] = 4.219309E-23; m[3126] = 8.227653E-24;
m[3128] = 3.291061E-24; m[3212] = 7.4E-20; m[3214] = 1.828367E-23; m[3216] = 5.485102E-24;
m[3218] = 3.656734E-24; m[3222] = 8.018E-11; m[3224] = 1.838582E-23; m[3226] = 5.485102E-24;
m[3228] = 3.656734E-24; m[3312] = 1.639E-10; m[3314] = 6.648608E-23; m[3322] = 2.9E-10;
m[3324] = 7.233101E-23; m[3334] = 8.21E-11; m[4112] = 2.991874E-22; m[4114] = 4.088274E-23;
m[4122] = 2.E-13; m[4132] = 1.12E-13; m[4212] = 3.999999E-22; m[4214] = 3.291061E-22;
m[4222] = 2.951624E-22; m[4224] = 4.417531E-23; m[4232] = 4.42E-13; m[4332] = 6.9E-14;
m[4412] = 3.335641E-13; m[4422] = 3.335641E-13; m[4432] = 3.335641E-13; m[5112] = 1.E-19;
m[5122] = 1.38E-12; m[5132] = 1.42E-12; m[5142] = 1.290893E-12; m[5212] = 1.E-19;
m[5222] = 1.E-19; m[5232] = 1.42E-12; m[5242] = 1.290893E-12; m[5312] = 1.E-19;
m[5322] = 1.E-19; m[5332] = 1.55E-12; m[5342] = 1.290893E-12; m[5442] = 1.290893E-12;
m[5512] = 1.290893E-12; m[5522] = 1.290893E-12; m[5532] = 1.290893E-12; m[5542] = 1.290893E-12;
m[10111] = 2.48382E-24; m[10113] = 4.635297E-24; m[10115] = 2.54136E-24; m[10211] = 2.48382E-24;
m[10213] = 4.635297E-24; m[10215] = 2.54136E-24; m[10223] = 1.828367E-24; m[10225] = 3.636531E-24;
m[10311] = 2.437823E-24; m[10313] = 7.313469E-24; m[10315] = 3.538775E-24;
m[10321] = 2.437823E-24; m[10323] = 7.313469E-24; m[10325] = 3.538775E-24;
m[10331] = 4.804469E-24; m[10411] = 4.38E-24; m[10413] = 3.29E-23; m[10421] = 4.38E-24;
m[10423] = 3.22653E-23; m[10431] = 6.5821E-22; m[10433] = 6.5821E-22; m[10441] = 6.453061E-23;
m[10511] = 4.39E-24; m[10513] = 1.65E-23; m[10521] = 4.39E-24; m[10523] = 1.65E-23;
m[10531] = 4.39E-24; m[10533] = 1.65E-23; m[11114] = 2.194041E-24; m[11116] = 1.828367E-24;
m[11212] = 1.880606E-24; m[11216] = 1.828367E-24; m[12112] = 2.194041E-24;
m[12114] = 2.194041E-24; m[12116] = 5.063171E-24; m[12126] = 1.828367E-24;
m[12212] = 2.194041E-24; m[12214] = 2.194041E-24; m[12216] = 5.063171E-24;
m[12224] = 2.194041E-24; m[12226] = 1.828367E-24; m[13112] = 6.582122E-24; m[13114] = 1.09702E-23;
m[13116] = 5.485102E-24; m[13122] = 1.316424E-23; m[13124] = 1.09702E-23; m[13126] = 6.928549E-24;
m[13212] = 6.582122E-24; m[13214] = 1.09702E-23; m[13216] = 5.485102E-24; m[13222] = 6.582122E-24;
m[13224] = 1.09702E-23; m[13226] = 5.485102E-24; m[13314] = 2.742551E-23;
m[13324] = 2.742551E-23; m[14122] = 1.828367E-22; m[20022] = 1.E+16; m[20113] = 1.567172E-24;
m[20213] = 1.567172E-24; m[20223] = 2.708692E-23; m[20313] = 3.782829E-24;
m[20315] = 2.384827E-24; m[20323] = 3.782829E-24; m[20325] = 2.384827E-24;
m[20333] = 1.198929E-23; m[20413] = 2.63E-24; m[20423] = 2.63E-24; m[20433] = 6.5821E-22;
m[20443] = 7.395643E-22; m[20513] = 2.63E-24; m[20523] = 2.63E-24; m[20533] = 2.63E-24;
m[21112] = 2.632849E-24; m[21114] = 3.291061E-24; m[21212] = 2.632849E-24;
m[21214] = 6.582122E-24; m[22112] = 4.388081E-24; m[22114] = 3.291061E-24;
m[22122] = 2.632849E-24; m[22124] = 6.582122E-24; m[22212] = 4.388081E-24;
m[22214] = 3.291061E-24; m[22222] = 2.632849E-24; m[22224] = 3.291061E-24;
m[23112] = 7.313469E-24; m[23114] = 2.991874E-24; m[23122] = 4.388081E-24;
m[23124] = 6.582122E-24; m[23126] = 3.291061E-24; m[23212] = 7.313469E-24;
m[23214] = 2.991874E-24; m[23222] = 7.313469E-24; m[23224] = 2.991874E-24;
m[30113] = 2.632849E-24; m[30213] = 2.632849E-24; m[30221] = 1.880606E-24;
m[30223] = 2.089563E-24; m[30313] = 2.056913E-24; m[30323] = 2.056913E-24;
m[30443] = 2.419898E-23; m[31114] = 1.880606E-24; m[31214] = 3.291061E-24;
m[32112] = 3.989164E-24; m[32114] = 1.880606E-24; m[32124] = 3.291061E-24;
m[32212] = 3.989164E-24; m[32214] = 1.880606E-24; m[32224] = 1.880606E-24;
m[33122] = 1.880606E-23; m[42112] = 6.582122E-24; m[42212] = 6.582122E-24;
m[43122] = 2.194041E-24; m[53122] = 4.388081E-24; m[100111] = 1.645531E-24;
m[100113] = 1.64553E-24; m[100211] = 1.645531E-24; m[100213] = 1.64553E-24;
m[100221] = 1.196749E-23; m[100223] = 3.061452E-24; m[100313] = 2.837122E-24;
m[100323] = 2.837122E-24; m[100331] = 4.459432E-25; m[100333] = 4.388081E-24;
m[100441] = 4.701516E-23; m[100443] = 2.076379E-21; m[100553] = 2.056913E-20;
m[200553] = 3.242425E-20; m[300553] = 3.210791E-23; m[9000111] = 8.776163E-24;
m[9000211] = 8.776163E-24; m[9000443] = 8.227652E-24; m[9000553] = 5.983747E-24;
m[9010111] = 3.164482E-24; m[9010211] = 3.164482E-24; m[9010221] = 9.403031E-24;
m[9010443] = 8.438618E-24; m[9010553] = 8.3318E-24; m[9020443] = 1.061633E-23;
m[9030221] = 6.038644E-24; m[9042413] = 2.07634E-21; m[9050225] = 1.394517E-24;
m[9060225] = 3.291061E-24; m[9080225] = 4.388081E-24; m[9090225] = 2.056913E-24;
m[9910445] = 2.07634E-21; m[9920443] = 2.07634E-21;
return true;
}
};
const int LHCB_2012_I1119400::_stablePDGIds[205] = {
311, 543, 545, 551, 555, 557, 1103, 2101, 2103, 2203, 3101, 3103, 3201, 3203, 3303,
4101, 4103, 4124, 4201, 4203, 4301, 4303, 4312, 4314, 4322, 4324, 4334, 4403, 4414,
4424, 4434, 4444, 5101, 5103, 5114, 5201, 5203, 5214, 5224, 5301, 5303, 5314, 5324,
5334, 5401, 5403, 5412, 5414, 5422, 5424, 5432, 5434, 5444, 5503, 5514, 5524, 5534,
5544, 5554, 10022, 10333, 10335, 10443, 10541, 10543, 10551, 10553, 10555, 11112,
12118, 12122, 12218, 12222, 13316, 13326, 20543, 20553, 20555, 23314, 23324, 30343,
30353, 30363, 30553, 33314, 33324, 41214, 42124, 52114, 52214, 100311, 100315, 100321,
100325, 100411, 100413, 100421, 100423, 100551, 100555, 100557, 110551, 110553, 110555,
120553, 120555, 130553, 200551, 200555, 210551, 210553, 220553, 1000001, 1000002,
1000003, 1000004, 1000005, 1000006, 1000011, 1000012, 1000013, 1000014, 1000015,
1000016, 1000021, 1000022, 1000023, 1000024, 1000025, 1000035, 1000037, 1000039,
2000001, 2000002, 2000003, 2000004, 2000005, 2000006, 2000011, 2000012, 2000013,
2000014, 2000015, 2000016, 3000111, 3000113, 3000211, 3000213, 3000221, 3000223,
3000331, 3100021, 3100111, 3100113, 3200111, 3200113, 3300113, 3400113, 4000001,
4000002, 4000011, 4000012, 5000039, 9000221, 9900012, 9900014, 9900016, 9900023,
9900024, 9900041, 9900042 };
// Plugin hook
DECLARE_RIVET_PLUGIN(LHCB_2012_I1119400);
}
diff --git a/analyses/pluginLHCb/LHCB_2012_I1208102.cc b/analyses/pluginLHCb/LHCB_2012_I1208102.cc.segfault
rename from analyses/pluginLHCb/LHCB_2012_I1208102.cc
rename to analyses/pluginLHCb/LHCB_2012_I1208102.cc.segfault
--- a/analyses/pluginLHCb/LHCB_2012_I1208102.cc
+++ b/analyses/pluginLHCb/LHCB_2012_I1208102.cc.segfault
@@ -1,79 +1,82 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/ZFinder.hh"
namespace Rivet {
/// Differential cross-sections of $\mathrm{Z}/\gamma^* \to e^{+}e^{-}$ vs rapidity and $\phi^*$
class LHCB_2012_I1208102 : public Analysis {
public:
/// Constructor
LHCB_2012_I1208102()
: Analysis("LHCB_2012_I1208102")
{ }
/// @name Analysis methods
//@{
/// Book histograms
void init() {
ZFinder zeefinder(FinalState(), Cuts::etaIn(2.0, 4.5) && Cuts::pT > 20*GeV, PID::ELECTRON, 60*GeV, 120*GeV);
declare(zeefinder, "ZeeFinder");
book(_h_sigma_vs_y ,2, 1, 1);
book(_h_sigma_vs_phi ,3, 1, 1);
}
/// Do the analysis
void analyze(const Event& e) {
const ZFinder& zeefinder = apply<ZFinder>(e, "ZeeFinder");
if (zeefinder.empty()) vetoEvent;
if (zeefinder.bosons().size() > 1)
MSG_WARNING("Found multiple (" << zeefinder.bosons().size() << ") Z -> e+ e- decays!");
// Z momenta
const FourMomentum& zee = zeefinder.bosons()[0].momentum();
+
+ if (zeefinder.constituents().size() < 2) vetoEvent;
+
const Particle& pozitron = zeefinder.constituents()[0];
const Particle& electron = zeefinder.constituents()[1];
// Calculation of the angular variable
const double diffphi = deltaPhi(pozitron, electron);
const double diffpsd = deltaEta(pozitron, electron);
const double accphi = M_PI - diffphi;
const double angular = tan(accphi/2) / cosh(diffpsd/2);
// Fill histograms
_h_sigma_vs_y->fill(zee.rapidity());
_h_sigma_vs_phi->fill(angular);
}
/// Finalize
void finalize() {
const double xs = crossSection()/picobarn;
scale(_h_sigma_vs_y, xs/sumOfWeights());
scale(_h_sigma_vs_phi, xs/sumOfWeights());
}
//@}
private:
/// @name Histograms
//@{
Histo1DPtr _h_sigma_vs_y, _h_sigma_vs_phi;
//@}
};
// The hook for the plugin system
DECLARE_RIVET_PLUGIN(LHCB_2012_I1208102);
}
diff --git a/analyses/pluginLHCb/LHCB_2013_I1208105.cc b/analyses/pluginLHCb/LHCB_2013_I1208105.cc
--- a/analyses/pluginLHCb/LHCB_2013_I1208105.cc
+++ b/analyses/pluginLHCb/LHCB_2013_I1208105.cc
@@ -1,231 +1,235 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/FinalState.hh"
#include "Rivet/Projections/ChargedFinalState.hh"
namespace Rivet {
class LHCB_2013_I1208105 : public Analysis {
public:
LHCB_2013_I1208105()
: Analysis("LHCB_2013_I1208105")
{ }
void init() {
// Projections
declare(FinalState(1.9, 4.9), "forwardFS");
declare(FinalState(-3.5,-1.5), "backwardFS");
declare(ChargedFinalState(1.9, 4.9), "forwardCFS");
declare(ChargedFinalState(-3.5,-1.5), "backwardCFS");
// Histos
book(_s_chEF_minbias, 1, 1, 1, true);
book(_s_chEF_hard, 2, 1, 1, true);
book(_s_chEF_diff, 3, 1, 1, true);
book(_s_chEF_nondiff, 4, 1, 1, true);
book(_s_totEF_minbias, 5, 1, 1, true);
book(_s_totEF_hard, 6, 1, 1, true);
book(_s_totEF_diff, 7, 1, 1, true);
book(_s_totEF_nondiff, 8, 1, 1, true);
// Temporary profiles and histos
/// @todo Convert to declared/registered temp histos
- _tp_chEF_minbias.reset(new YODA::Profile1D(refData(1,1,1)));
- _tp_chEF_hard.reset(new YODA::Profile1D(refData(2,1,1)));
- _tp_chEF_diff.reset(new YODA::Profile1D(refData(3,1,1)));
- _tp_chEF_nondiff.reset(new YODA::Profile1D(refData(4,1,1)));
- _tp_totEF_minbias.reset(new YODA::Profile1D(refData(5,1,1)));
- _tp_totEF_hard.reset(new YODA::Profile1D(refData(6,1,1)));
- _tp_totEF_diff.reset(new YODA::Profile1D(refData(7,1,1)));
- _tp_totEF_nondiff.reset(new YODA::Profile1D(refData(8,1,1)));
- //
- _th_chN_minbias.reset(new YODA::Histo1D(refData(1,1,1)));
- _th_chN_hard.reset(new YODA::Histo1D(refData(2,1,1)));
- _th_chN_diff.reset(new YODA::Histo1D(refData(3,1,1)));
- _th_chN_nondiff.reset(new YODA::Histo1D(refData(4,1,1)));
- _th_totN_minbias.reset(new YODA::Histo1D(refData(5,1,1)));
- _th_totN_hard.reset(new YODA::Histo1D(refData(6,1,1)));
- _th_totN_diff.reset(new YODA::Histo1D(refData(7,1,1)));
- _th_totN_nondiff.reset(new YODA::Histo1D(refData(8,1,1)));
+ book(_tp_chEF_minbias, "TMP/chEF_minbias", refData(1,1,1));
+ book(_tp_chEF_hard, "TMP/chEF_hard", refData(2,1,1));
+ book(_tp_chEF_diff, "TMP/chEF_diff", refData(3,1,1));
+ book(_tp_chEF_nondiff, "TMP/chEF_nondiff", refData(4,1,1));
+ book(_tp_totEF_minbias, "TMP/totEF_minbias", refData(5,1,1));
+ book(_tp_totEF_hard, "TMP/totEF_hard", refData(6,1,1));
+ book(_tp_totEF_diff, "TMP/totEF_diff", refData(7,1,1));
+ book(_tp_totEF_nondiff, "TMP/totEF_nondiff", refData(8,1,1));
+
+ book(_th_chN_minbias, "TMP/chN_minbias", refData(1,1,1));
+ book(_th_chN_hard, "TMP/chN_hard", refData(2,1,1));
+ book(_th_chN_diff, "TMP/chN_diff", refData(3,1,1));
+ book(_th_chN_nondiff, "TMP/chN_nondiff", refData(4,1,1));
+ book(_th_totN_minbias, "TMP/totN_minbias", refData(5,1,1));
+ book(_th_totN_hard, "TMP/totN_hard", refData(6,1,1));
+ book(_th_totN_diff, "TMP/totN_diff", refData(7,1,1));
+ book(_th_totN_nondiff, "TMP/totN_nondiff", refData(8,1,1));
// Counters
- _mbSumW = 0.0; _hdSumW = 0.0; _dfSumW = 0.0; _ndSumW = 0.0;
- _mbchSumW = 0.0; _hdchSumW = 0.0; _dfchSumW = 0.0; _ndchSumW = 0.0;
+ book(_mbSumW, "TMP/mbSumW");
+ book(_hdSumW, "TMP/hdSumW");
+ book(_dfSumW, "TMP/dfSumW");
+ book(_ndSumW, "TMP/ndSumW");
+ book(_mbchSumW, "TMP/mbchSumW");
+ book(_hdchSumW, "TMP/hdchSumW");
+ book(_dfchSumW, "TMP/dfchSumW");
+ book(_ndchSumW, "TMP/ndchSumW");
}
/// Perform the per-event analysis
void analyze(const Event& event) {
- const double weight = 1.0;
-
const FinalState& ffs = apply<FinalState>(event, "forwardFS");
const FinalState& bfs = apply<FinalState>(event, "backwardFS");
const ChargedFinalState& fcfs = apply<ChargedFinalState>(event, "forwardCFS");
const ChargedFinalState& bcfs = apply<ChargedFinalState>(event, "backwardCFS");
// Veto this event completely if there are no forward *charged* particles
if (fcfs.empty()) vetoEvent;
// Charged and neutral version
{
// Decide empirically if this is a "hard" or "diffractive" event
bool ishardEvt = false;
foreach (const Particle& p, ffs.particles()) {
if (p.pT() > 3.0*GeV) { ishardEvt = true; break; }
}
// Decide empirically if this is a "diffractive" event
/// @todo Can be "diffractive" *and* "hard"?
bool isdiffEvt = (bfs.size() == 0);
// Update event-type weight counters
- _mbSumW += weight;
- (isdiffEvt ? _dfSumW : _ndSumW) += weight;
- if (ishardEvt) _hdSumW += weight;
+ _mbSumW->fill();
+ (isdiffEvt ? _dfSumW : _ndSumW)->fill();
+ if (ishardEvt) _hdSumW->fill();
// Plot energy flow
foreach (const Particle& p, ffs.particles()) {
const double eta = p.eta();
const double energy = p.E();
- _tp_totEF_minbias->fill(eta, energy, weight);
- _th_totN_minbias->fill(eta, weight);
+ _tp_totEF_minbias->fill(eta, energy);
+ _th_totN_minbias->fill(eta);
if (ishardEvt) {
- _tp_totEF_hard->fill(eta, energy, weight);
- _th_totN_hard->fill(eta, weight);
+ _tp_totEF_hard->fill(eta, energy);
+ _th_totN_hard->fill(eta);
}
if (isdiffEvt) {
- _tp_totEF_diff->fill(eta, energy, weight);
- _th_totN_diff->fill(eta, weight);
+ _tp_totEF_diff->fill(eta, energy);
+ _th_totN_diff->fill(eta);
} else {
- _tp_totEF_nondiff->fill(eta, energy, weight);
- _th_totN_nondiff->fill(eta, weight);
+ _tp_totEF_nondiff->fill(eta, energy);
+ _th_totN_nondiff->fill(eta);
}
}
}
// Charged-only version
{
bool ishardEvt = false;
foreach (const Particle& p, fcfs.particles()) {
if (p.pT() > 3.0*GeV) { ishardEvt = true; break; }
}
// Decide empirically if this is a "diffractive" event
/// @todo Can be "diffractive" *and* "hard"?
bool isdiffEvt = (bcfs.size() == 0);
// Update event-type weight counters
- _mbchSumW += weight;
- (isdiffEvt ? _dfchSumW : _ndchSumW) += weight;
- if (ishardEvt) _hdchSumW += weight;
+ _mbchSumW->fill();
+ (isdiffEvt ? _dfchSumW : _ndchSumW)->fill();
+ if (ishardEvt) _hdchSumW->fill();
// Plot energy flow
foreach (const Particle& p, fcfs.particles()) {
const double eta = p.eta();
const double energy = p.E();
- _tp_chEF_minbias->fill(eta, energy, weight);
- _th_chN_minbias->fill(eta, weight);
+ _tp_chEF_minbias->fill(eta, energy);
+ _th_chN_minbias->fill(eta);
if (ishardEvt) {
- _tp_chEF_hard->fill(eta, energy, weight);
- _th_chN_hard->fill(eta, weight);
+ _tp_chEF_hard->fill(eta, energy);
+ _th_chN_hard->fill(eta);
}
if (isdiffEvt) {
- _tp_chEF_diff->fill(eta, energy, weight);
- _th_chN_diff->fill(eta, weight);
+ _tp_chEF_diff->fill(eta, energy);
+ _th_chN_diff->fill(eta);
} else {
- _tp_chEF_nondiff->fill(eta, energy, weight);
- _th_chN_nondiff->fill(eta, weight);
+ _tp_chEF_nondiff->fill(eta, energy);
+ _th_chN_nondiff->fill(eta);
}
}
}
}
void finalize() {
for (size_t i = 0; i < _s_totEF_minbias->numPoints(); ++i) {
const double val = _tp_totEF_minbias->bin(i).mean() * _th_totN_minbias->bin(i).height();
const double err = (_tp_totEF_minbias->bin(i).mean() * _th_totN_minbias->bin(i).heightErr() +
_tp_totEF_minbias->bin(i).stdErr() * _th_totN_minbias->bin(i).height());
_s_totEF_minbias->point(i).setY(val/_mbSumW, err/_mbSumW);
}
for (size_t i = 0; i < _s_totEF_hard->numPoints(); ++i) {
const double val = _tp_totEF_hard->bin(i).mean() * _th_totN_hard->bin(i).height();
const double err = (_tp_totEF_hard->bin(i).mean() * _th_totN_hard->bin(i).heightErr() +
_tp_totEF_hard->bin(i).stdErr() * _th_totN_hard->bin(i).height());
_s_totEF_hard->point(i).setY(val/_hdSumW, err/_hdSumW);
}
for (size_t i = 0; i < _s_totEF_diff->numPoints(); ++i) {
const double val = _tp_totEF_diff->bin(i).mean() * _th_totN_diff->bin(i).height();
const double err = (_tp_totEF_diff->bin(i).mean() * _th_totN_diff->bin(i).heightErr() +
_tp_totEF_diff->bin(i).stdErr() * _th_totN_diff->bin(i).height());
_s_totEF_diff->point(i).setY(val/_dfSumW, err/_dfSumW);
}
for (size_t i = 0; i < _s_totEF_nondiff->numPoints(); ++i) {
const double val = _tp_totEF_nondiff->bin(i).mean() * _th_totN_nondiff->bin(i).height();
const double err = (_tp_totEF_nondiff->bin(i).mean() * _th_totN_nondiff->bin(i).heightErr() +
_tp_totEF_nondiff->bin(i).stdErr() * _th_totN_nondiff->bin(i).height());
_s_totEF_nondiff->point(i).setY(val/_ndSumW, err/_ndSumW);
}
for (size_t i = 0; i < _s_chEF_minbias->numPoints(); ++i) {
const double val = _tp_chEF_minbias->bin(i).mean() * _th_chN_minbias->bin(i).height();
const double err = (_tp_chEF_minbias->bin(i).mean() * _th_chN_minbias->bin(i).heightErr() +
_tp_chEF_minbias->bin(i).stdErr() * _th_chN_minbias->bin(i).height());
_s_chEF_minbias->point(i).setY(val/_mbchSumW, err/_mbchSumW);
}
for (size_t i = 0; i < _s_chEF_hard->numPoints(); ++i) {
const double val = _tp_chEF_hard->bin(i).mean() * _th_chN_hard->bin(i).height();
const double err = (_tp_chEF_hard->bin(i).mean() * _th_chN_hard->bin(i).heightErr() +
_tp_chEF_hard->bin(i).stdErr() * _th_chN_hard->bin(i).height());
_s_chEF_hard->point(i).setY(val/_hdchSumW, err/_hdchSumW);
}
for (size_t i = 0; i < _s_chEF_diff->numPoints(); ++i) {
const double val = _tp_chEF_diff->bin(i).mean() * _th_chN_diff->bin(i).height();
const double err = (_tp_chEF_diff->bin(i).mean() * _th_chN_diff->bin(i).heightErr() +
_tp_chEF_diff->bin(i).stdErr() * _th_chN_diff->bin(i).height());
_s_chEF_diff->point(i).setY(val/_dfchSumW, err/_dfchSumW);
}
for (size_t i = 0; i < _s_chEF_nondiff->numPoints(); ++i) {
const double val = _tp_chEF_nondiff->bin(i).mean() * _th_chN_nondiff->bin(i).height();
const double err = (_tp_chEF_nondiff->bin(i).mean() * _th_chN_nondiff->bin(i).heightErr() +
_tp_chEF_nondiff->bin(i).stdErr() * _th_chN_nondiff->bin(i).height());
_s_chEF_nondiff->point(i).setY(val/_ndchSumW, err/_ndchSumW);
}
}
private:
/// @name Histograms and counters
///
/// @note Histograms correspond to charged and total EF for each class of events:
/// minimum bias, hard scattering, diffractive enriched and non-diffractive enriched.
//@{
// Scatters to be filled in finalize with 1/d_eta <N(eta)><E(eta)>
Scatter2DPtr _s_totEF_minbias, _s_totEF_hard, _s_totEF_diff, _s_totEF_nondiff;
Scatter2DPtr _s_chEF_minbias, _s_chEF_hard, _s_chEF_diff, _s_chEF_nondiff;
// Temp profiles containing <E(eta)>
- std::shared_ptr<YODA::Profile1D> _tp_totEF_minbias, _tp_totEF_hard, _tp_totEF_diff, _tp_totEF_nondiff;
- std::shared_ptr<YODA::Profile1D> _tp_chEF_minbias, _tp_chEF_hard, _tp_chEF_diff, _tp_chEF_nondiff;
+ Profile1DPtr _tp_totEF_minbias, _tp_totEF_hard, _tp_totEF_diff, _tp_totEF_nondiff;
+ Profile1DPtr _tp_chEF_minbias, _tp_chEF_hard, _tp_chEF_diff, _tp_chEF_nondiff;
// Temp profiles containing <N(eta)>
- std::shared_ptr<YODA::Histo1D> _th_totN_minbias, _th_totN_hard, _th_totN_diff, _th_totN_nondiff;
- std::shared_ptr<YODA::Histo1D> _th_chN_minbias, _th_chN_hard, _th_chN_diff, _th_chN_nondiff;
+ Histo1DPtr _th_totN_minbias, _th_totN_hard, _th_totN_diff, _th_totN_nondiff;
+ Histo1DPtr _th_chN_minbias, _th_chN_hard, _th_chN_diff, _th_chN_nondiff;
// Sums of weights (~ #events) in each event class
- double _mbSumW, _hdSumW, _dfSumW, _ndSumW;
- double _mbchSumW, _hdchSumW, _dfchSumW, _ndchSumW;
+ CounterPtr _mbSumW, _hdSumW, _dfSumW, _ndSumW;
+ CounterPtr _mbchSumW, _hdchSumW, _dfchSumW, _ndchSumW;
//@}
};
// Hook for the plugin system
DECLARE_RIVET_PLUGIN(LHCB_2013_I1208105);
}
diff --git a/analyses/pluginLHCb/LHCB_2014_I1281685.cc b/analyses/pluginLHCb/LHCB_2014_I1281685.cc
--- a/analyses/pluginLHCb/LHCB_2014_I1281685.cc
+++ b/analyses/pluginLHCb/LHCB_2014_I1281685.cc
@@ -1,1178 +1,1177 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/ChargedFinalState.hh"
namespace Rivet {
/// Charged particle multiplicities and densities in $pp$ collisions at $\sqrt{s} = 7$ TeV
class LHCB_2014_I1281685 : public Analysis {
public:
/// @name Constructors etc.
//@{
/// Constructor
LHCB_2014_I1281685()
: Analysis("LHCB_2014_I1281685"),
_p_min(2.0),
_pt_min(0.2),
_eta_min(2.0),
_eta_max(4.8),
_maxlft(1.0e-11)
{ }
//@}
/// @name Analysis methods
//@{
/// Book histograms and initialise projections before the run
void init() {
fillMap(_partLftMap);
// Projections
declare(ChargedFinalState(_eta_min, _eta_max, _pt_min*GeV), "CFS");
// Book histograms
book(_h_mult_total ,"d03-x01-y01", 50, 0.5, 50.5);
book(_h_mult_eta[0] ,"d04-x01-y01", 21, -0.5, 20.5); //eta=[2.0,2.5]
book(_h_mult_eta[1] ,"d04-x01-y02", 21, -0.5, 20.5); //eta=[2.5,3.0]
book(_h_mult_eta[2] ,"d04-x01-y03", 21, -0.5, 20.5); //eta=[3.0,3.5]
book(_h_mult_eta[3] ,"d04-x01-y04", 21, -0.5, 20.5); //eta=[3.5,4.0]
book(_h_mult_eta[4] ,"d04-x01-y05", 21, -0.5, 20.5); //eta=[4.0,4.5]
book(_h_mult_pt[0] ,"d05-x01-y01", 21, -0.5, 20.5); //pT=[0.2,0.3]GeV
book(_h_mult_pt[1] ,"d05-x01-y02", 21, -0.5, 20.5); //pT=[0.3,0.4]GeV
book(_h_mult_pt[2] ,"d05-x01-y03", 21, -0.5, 20.5); //pT=[0.4,0.6]GeV
book(_h_mult_pt[3] ,"d05-x01-y04", 21, -0.5, 20.5); //pT=[0.6,1.0]GeV
book(_h_mult_pt[4] ,"d05-x01-y05", 21, -0.5, 20.5); //pT=[1.0,2.0]GeV
book(_h_dndeta ,"d01-x01-y01", 14, 2.0, 4.8); //eta=[2,4.8]
book(_h_dndpt ,"d02-x01-y01", 18, 0.2, 2.0); //pT =[0,2]GeV
// Counters
- _sumW = 0;
+ book(_sumW, "TMP/sumW");
}
/// Perform the per-event analysis
void analyze(const Event& event) {
// Variable to store multiplicities per event
int LHCbcountAll = 0; //count particles fulfiling all requirements
int LHCbcountEta[8] = {0,0,0,0,0,0,0,0}; //count per eta-bin
int LHCbcountPt[7] = {0,0,0,0,0,0,0}; //count per pT-bin
vector<double> val_dNdEta;
vector<double> val_dNdPt;
val_dNdEta.clear();
val_dNdPt.clear();
const ChargedFinalState& cfs = apply<ChargedFinalState>(event, "CFS");
foreach (const Particle& p, cfs.particles()) {
int id = p.pdgId();
// continue if particle is not a pion, kaon, proton, muon or electron
if ( !( (abs(id) == 211) || (abs(id) == 321) || (abs(id) == 2212) || (abs(id) == 13) || (abs(id) == 11)) ) {
continue;
}
const FourMomentum& qmom = p.momentum();
const double eta = p.momentum().eta();
const double pT = p.momentum().pT();
//minimum momentum
if (qmom.p3().mod() < _p_min) continue;
//minimum tr. momentum
if (pT < _pt_min) continue;
//eta range
if ((eta < _eta_min) || (eta > _eta_max)) continue;
/* Select only prompt particles via lifetime */
//Sum of all mother lifetimes (PDG lifetime) < 10ps
double ancestors_sumlft = getAncestorSumLifetime(p);
if( (ancestors_sumlft > _maxlft) || (ancestors_sumlft < 0) ) continue;
//after all cuts;
LHCbcountAll++; //count particles in whole kin. range
//in eta bins
if( eta >2.0 && eta <= 2.5) LHCbcountEta[0]++;
if( eta >2.5 && eta <= 3.0) LHCbcountEta[1]++;
if( eta >3.0 && eta <= 3.5) LHCbcountEta[2]++;
if( eta >3.5 && eta <= 4.0) LHCbcountEta[3]++;
if( eta >4.0 && eta <= 4.5) LHCbcountEta[4]++;
if( eta >2.0 && eta <= 4.8) LHCbcountEta[5]++; //cross-check
//in pT bins
if( pT > 0.2 && pT <= 0.3) LHCbcountPt[0]++;
if( pT > 0.3 && pT <= 0.4) LHCbcountPt[1]++;
if( pT > 0.4 && pT <= 0.6) LHCbcountPt[2]++;
if( pT > 0.6 && pT <= 1.0) LHCbcountPt[3]++;
if( pT > 1.0 && pT <= 2.0) LHCbcountPt[4]++;
if( pT > 0.2) LHCbcountPt[5]++; //cross-check
//particle densities -> need proper normalization (finalize)
val_dNdPt.push_back( pT );
val_dNdEta.push_back( eta );
}//end foreach
// Fill histograms only, if at least 1 particle pre event was within the
// kinematic range of the analysis!
if (LHCbcountAll) {
- const double weight = 1.0;
- _sumW += weight;
+ _sumW->fill();
- _h_mult_total->fill(LHCbcountAll, weight);
+ _h_mult_total->fill(LHCbcountAll);
- _h_mult_eta[0]->fill(LHCbcountEta[0], weight);
- _h_mult_eta[1]->fill(LHCbcountEta[1], weight);
- _h_mult_eta[2]->fill(LHCbcountEta[2], weight);
- _h_mult_eta[3]->fill(LHCbcountEta[3], weight);
- _h_mult_eta[4]->fill(LHCbcountEta[4], weight);
+ _h_mult_eta[0]->fill(LHCbcountEta[0]);
+ _h_mult_eta[1]->fill(LHCbcountEta[1]);
+ _h_mult_eta[2]->fill(LHCbcountEta[2]);
+ _h_mult_eta[3]->fill(LHCbcountEta[3]);
+ _h_mult_eta[4]->fill(LHCbcountEta[4]);
- _h_mult_pt[0]->fill(LHCbcountPt[0], weight);
- _h_mult_pt[1]->fill(LHCbcountPt[1], weight);
- _h_mult_pt[2]->fill(LHCbcountPt[2], weight);
- _h_mult_pt[3]->fill(LHCbcountPt[3], weight);
- _h_mult_pt[4]->fill(LHCbcountPt[4], weight);
+ _h_mult_pt[0]->fill(LHCbcountPt[0]);
+ _h_mult_pt[1]->fill(LHCbcountPt[1]);
+ _h_mult_pt[2]->fill(LHCbcountPt[2]);
+ _h_mult_pt[3]->fill(LHCbcountPt[3]);
+ _h_mult_pt[4]->fill(LHCbcountPt[4]);
for (size_t part = 0; part < val_dNdEta.size(); part++)
- _h_dndeta->fill(val_dNdEta[part], weight);
+ _h_dndeta->fill(val_dNdEta[part]);
for (size_t part = 0; part < val_dNdPt.size(); part++)
- _h_dndpt->fill(val_dNdPt[part], weight);
+ _h_dndpt->fill(val_dNdPt[part]);
}
}
/// Normalise histograms etc., after the run
void finalize() {
const double scalefactor = 1.0/_sumW; // normalize multiplicity histograms by nEvents
const double scale1k = 1000.; // to match '10^3' scale in reference histograms
scale( _h_dndeta, scalefactor );
scale( _h_dndpt, scalefactor*0.1 ); //additional factor 0.1 for [0.1 GeV/c]
scale( _h_mult_total, scalefactor*scale1k);
_h_mult_eta[0]->scaleW( scalefactor*scale1k );
_h_mult_eta[1]->scaleW( scalefactor*scale1k );
_h_mult_eta[2]->scaleW( scalefactor*scale1k );
_h_mult_eta[3]->scaleW( scalefactor*scale1k );
_h_mult_eta[4]->scaleW( scalefactor*scale1k );
_h_mult_pt[0]->scaleW( scalefactor*scale1k );
_h_mult_pt[1]->scaleW( scalefactor*scale1k );
_h_mult_pt[2]->scaleW( scalefactor*scale1k );
_h_mult_pt[3]->scaleW( scalefactor*scale1k );
_h_mult_pt[4]->scaleW( scalefactor*scale1k );
}
//@}
private:
// Get mean PDG lifetime for particle with PID
double getLifetime(int pid) {
double lft = 0.;
map<int, double>::iterator pPartLft = _partLftMap.find(pid);
if (pPartLft != _partLftMap.end()) {
lft = (*pPartLft).second;
} else {
// allow identifying missing life times only in debug mode
MSG_DEBUG("Could not determine lifetime for particle with PID " << pid << "... Assume non-prompt particle");
lft = -1;
}
return lft;
}
// Get sum of all ancestor particles
const double getAncestorSumLifetime(const Particle& p) {
double lftSum = 0.;
double plft = 0.;
const GenParticle* part = p.genParticle();
if ( 0 == part ) return -1;
const GenVertex* ivtx = part->production_vertex();
while(ivtx) {
if (ivtx->particles_in_size() < 1) { lftSum = -1.; break; };
const GenVertex::particles_in_const_iterator iPart_invtx = ivtx->particles_in_const_begin();
part = (*iPart_invtx);
if ( !(part) ) { lftSum = -1.; break; };
ivtx = part->production_vertex();
if ( (part->pdg_id() == 2212) || !(ivtx) ) break; // reached beam
plft = getLifetime(part->pdg_id());
if (plft < 0.) { lftSum = -1.; break; };
lftSum += plft;
}
return (lftSum);
}
/// Hard-coded map linking PDG ID with PDG lifetime[s] (converted from ParticleTable.txt)
bool fillMap(map<int, double>& m) {
// PDGID = LIFETIME
m[22] = 1.000000e+016;
m[-11] = 1.000000e+016;
m[11] = 1.000000e+016;
m[12] = 1.000000e+016;
m[-13] = 2.197036e-006;
m[13] = 2.197036e-006;
m[111] = 8.438618e-017;
m[211] = 2.603276e-008;
m[-211] = 2.603276e-008;
m[130] = 5.174624e-008;
m[321] = 1.238405e-008;
m[-321] = 1.238405e-008;
m[2112] = 885.646128;
m[2212] = 1.000000e+016;
m[-2212] = 1.000000e+016;
m[310] = 8.934603e-011;
m[221] = 5.578070e-019;
m[3122] = 2.631796e-010;
m[3222] = 8.018178e-011;
m[3212] = 7.395643e-020;
m[3112] = 1.479129e-010;
m[3322] = 2.899613e-010;
m[3312] = 1.637344e-010;
m[3334] = 8.207135e-011;
m[-2112] = 885.646128;
m[-3122] = 2.631796e-010;
m[-3222] = 8.018178e-011;
m[-3212] = 7.395643e-020;
m[-3112] = 1.479129e-010;
m[-3322] = 2.899613e-010;
m[-3312] = 1.637344e-010;
m[-3334] = 8.207135e-011;
m[113] = 4.411610e-024;
m[213] = 4.411610e-024;
m[-213] = 4.411610e-024;
m[223] = 7.798723e-023;
m[333] = 1.545099e-022;
m[323] = 1.295693e-023;
m[-323] = 1.295693e-023;
m[313] = 1.298249e-023;
m[-313] = 1.298249e-023;
m[20213] = 1.500000e-024;
m[-20213] = 1.500000e-024;
m[450000000] = 1.000000e+015;
m[460000000] = 1.000000e+015;
m[470000000] = 1.000000e+015;
m[480000000] = 1.000000e+015;
m[490000000] = 1.000000e+015;
m[20022] = 1.000000e+016;
m[-15] = 2.906014e-013;
m[15] = 2.906014e-013;
m[24] = 3.104775e-025;
m[-24] = 3.104775e-025;
m[23] = 2.637914e-025;
m[411] = 1.051457e-012;
m[-411] = 1.051457e-012;
m[421] = 4.116399e-013;
m[-421] = 4.116399e-013;
m[431] = 4.904711e-013;
m[-431] = 4.904711e-013;
m[4122] = 1.994582e-013;
m[-4122] = 1.994582e-013;
m[443] = 7.565657e-021;
m[413] = 6.856377e-021;
m[-413] = 6.856377e-021;
m[423] = 1.000003e-019;
m[-423] = 1.000003e-019;
m[433] = 1.000003e-019;
m[-433] = 1.000003e-019;
m[521] = 1.671000e-012;
m[-521] = 1.671000e-012;
m[511] = 1.536000e-012;
m[-511] = 1.536000e-012;
m[531] = 1.461000e-012;
m[-531] = 1.461000e-012;
m[541] = 4.600000e-013;
m[-541] = 4.600000e-013;
m[5122] = 1.229000e-012;
m[-5122] = 1.229000e-012;
m[4112] = 4.388081e-022;
m[-4112] = 4.388081e-022;
m[4212] = 3.999999e-022;
m[-4212] = 3.999999e-022;
m[4222] = 3.291060e-022;
m[-4222] = 3.291060e-022;
m[25] = 9.400000e-026;
m[35] = 9.400000e-026;
m[36] = 9.400000e-026;
m[37] = 9.400000e-026;
m[-37] = 9.400000e-026;
m[4312] = 9.800002e-014;
m[-4312] = 9.800002e-014;
m[4322] = 3.500001e-013;
m[-4322] = 3.500001e-013;
m[4332] = 6.453061e-014;
m[-4332] = 6.453061e-014;
m[4132] = 9.824063e-014;
m[-4132] = 9.824063e-014;
m[4232] = 4.417532e-013;
m[-4232] = 4.417532e-013;
m[5222] = 1.000000e-019;
m[-5222] = 1.000000e-019;
m[5212] = 1.000000e-019;
m[-5212] = 1.000000e-019;
m[5112] = 1.000000e-019;
m[-5112] = 1.000000e-019;
m[5312] = 1.000000e-019;
m[-5312] = 1.000000e-019;
m[5322] = 1.000000e-019;
m[-5322] = 1.000000e-019;
m[5332] = 1.550000e-012;
m[-5332] = 1.550000e-012;
m[5132] = 1.390000e-012;
m[-5132] = 1.390000e-012;
m[5232] = 1.390000e-012;
m[-5232] = 1.390000e-012;
m[100443] = 2.194041e-021;
m[331] = 3.258476e-021;
m[441] = 4.113826e-023;
m[10441] = 4.063038e-023;
m[20443] = 7.154480e-022;
m[445] = 3.164482e-022;
m[9000111] = 1.149997e-023;
m[9000211] = 1.149997e-023;
m[-9000211] = 1.149997e-023;
m[20113] = 1.500000e-024;
m[115] = 6.151516e-024;
m[215] = 6.151516e-024;
m[-215] = 6.151516e-024;
m[10323] = 7.313469e-024;
m[-10323] = 7.313469e-024;
m[10313] = 7.313469e-024;
m[-10313] = 7.313469e-024;
m[20323] = 3.782829e-024;
m[-20323] = 3.782829e-024;
m[20313] = 3.782829e-024;
m[-20313] = 3.782829e-024;
m[10321] = 2.238817e-024;
m[-10321] = 2.238817e-024;
m[10311] = 2.238817e-024;
m[-10311] = 2.238817e-024;
m[325] = 6.682357e-024;
m[-325] = 6.682357e-024;
m[315] = 6.038644e-024;
m[-315] = 6.038644e-024;
m[10411] = 4.380000e-024;
m[20413] = 2.630000e-024;
m[10413] = 3.290000e-023;
m[-415] = 2.632849e-023;
m[-10411] = 4.380000e-024;
m[-20413] = 2.630000e-024;
m[-10413] = 3.290000e-023;
m[415] = 2.632849e-023;
m[10421] = 4.380000e-024;
m[20423] = 2.630000e-024;
m[10423] = 3.482604e-023;
m[-425] = 2.861792e-023;
m[-10421] = 4.380000e-024;
m[-20423] = 2.630000e-024;
m[-10423] = 3.482604e-023;
m[425] = 2.861792e-023;
m[10431] = 6.582100e-022;
m[20433] = 6.582100e-022;
m[10433] = 6.582100e-022;
m[435] = 4.388100e-023;
m[-10431] = 6.582100e-022;
m[-20433] = 6.582100e-022;
m[-10433] = 6.582100e-022;
m[-435] = 4.388100e-023;
m[2224] = 5.485102e-024;
m[2214] = 5.485102e-024;
m[2114] = 5.485102e-024;
m[1114] = 5.485102e-024;
m[-2224] = 5.485102e-024;
m[-2214] = 5.485102e-024;
m[-2114] = 5.485102e-024;
m[-1114] = 5.485102e-024;
m[-523] = 1.000019e-019;
m[523] = 1.000019e-019;
m[513] = 1.000019e-019;
m[-513] = 1.000019e-019;
m[533] = 1.000000e-019;
m[-533] = 1.000000e-019;
m[10521] = 4.390000e-024;
m[20523] = 2.630000e-024;
m[10523] = 1.650000e-023;
m[525] = 1.310000e-023;
m[-10521] = 4.390000e-024;
m[-20523] = 2.630000e-024;
m[-10523] = 1.650000e-023;
m[-525] = 1.310000e-023;
m[10511] = 4.390000e-024;
m[20513] = 2.630000e-024;
m[10513] = 1.650000e-023;
m[515] = 1.310000e-023;
m[-10511] = 4.390000e-024;
m[-20513] = 2.630000e-024;
m[-10513] = 1.650000e-023;
m[-515] = 1.310000e-023;
m[10531] = 4.390000e-024;
m[20533] = 2.630000e-024;
m[10533] = 1.650000e-023;
m[535] = 1.310000e-023;
m[-10531] = 4.390000e-024;
m[-20533] = 2.630000e-024;
m[-10533] = 1.650000e-023;
m[-535] = 1.310000e-023;
m[14] = 1.000000e+016;
m[-14] = 1.000000e+016;
m[-12] = 1.000000e+016;
m[1] = 0.000000e+000;
m[-1] = 0.000000e+000;
m[2] = 0.000000e+000;
m[-2] = 0.000000e+000;
m[3] = 0.000000e+000;
m[-3] = 0.000000e+000;
m[4] = 0.000000e+000;
m[-4] = 0.000000e+000;
m[5] = 0.000000e+000;
m[-5] = 0.000000e+000;
m[6] = 4.707703e-025;
m[-6] = 4.707703e-025;
m[7] = 0.000000e+000;
m[-7] = 0.000000e+000;
m[8] = 0.000000e+000;
m[-8] = 0.000000e+000;
m[16] = 1.000000e+016;
m[-16] = 1.000000e+016;
m[17] = 0.000000e+000;
m[-17] = 0.000000e+000;
m[18] = 0.000000e+000;
m[-18] = 0.000000e+000;
m[21] = 0.000000e+000;
m[32] = 0.000000e+000;
m[33] = 0.000000e+000;
m[34] = 0.000000e+000;
m[-34] = 0.000000e+000;
m[39] = 0.000000e+000;
m[41] = 0.000000e+000;
m[-41] = 0.000000e+000;
m[42] = 0.000000e+000;
m[-42] = 0.000000e+000;
m[43] = 0.000000e+000;
m[44] = 0.000000e+000;
m[-44] = 0.000000e+000;
m[81] = 0.000000e+000;
m[82] = 0.000000e+000;
m[-82] = 0.000000e+000;
m[83] = 0.000000e+000;
m[84] = 3.335641e-013;
m[-84] = 3.335641e-013;
m[85] = 1.290893e-012;
m[-85] = 1.290893e-012;
m[86] = 0.000000e+000;
m[-86] = 0.000000e+000;
m[87] = 0.000000e+000;
m[-87] = 0.000000e+000;
m[88] = 0.000000e+000;
m[90] = 0.000000e+000;
m[91] = 0.000000e+000;
m[92] = 0.000000e+000;
m[93] = 0.000000e+000;
m[94] = 0.000000e+000;
m[95] = 0.000000e+000;
m[96] = 0.000000e+000;
m[97] = 0.000000e+000;
m[98] = 0.000000e+000;
m[99] = 0.000000e+000;
m[117] = 4.088275e-024;
m[119] = 1.828367e-024;
m[217] = 4.088275e-024;
m[-217] = 4.088275e-024;
m[219] = 1.828367e-024;
m[-219] = 1.828367e-024;
m[225] = 3.555982e-024;
m[227] = 3.917930e-024;
m[229] = 3.392846e-024;
m[311] = 1.000000e+016;
m[-311] = 1.000000e+016;
m[317] = 4.139699e-024;
m[-317] = 4.139699e-024;
m[319] = 3.324304e-024;
m[-319] = 3.324304e-024;
m[327] = 4.139699e-024;
m[-327] = 4.139699e-024;
m[329] = 3.324304e-024;
m[-329] = 3.324304e-024;
m[335] = 8.660687e-024;
m[337] = 7.565657e-024;
m[543] = 0.000000e+000;
m[-543] = 0.000000e+000;
m[545] = 0.000000e+000;
m[-545] = 0.000000e+000;
m[551] = 0.000000e+000;
m[553] = 1.253738e-020;
m[555] = 1.000000e+016;
m[557] = 0.000000e+000;
m[-450000000] = 0.000000e+000;
m[-490000000] = 0.000000e+000;
m[-460000000] = 0.000000e+000;
m[-470000000] = 0.000000e+000;
m[1103] = 0.000000e+000;
m[-1103] = 0.000000e+000;
m[1112] = 4.388081e-024;
m[-1112] = 4.388081e-024;
m[1116] = 1.880606e-024;
m[-1116] = 1.880606e-024;
m[1118] = 2.194041e-024;
m[-1118] = 2.194041e-024;
m[1212] = 4.388081e-024;
m[-1212] = 4.388081e-024;
m[1214] = 5.485102e-024;
m[-1214] = 5.485102e-024;
m[1216] = 1.880606e-024;
m[-1216] = 1.880606e-024;
m[1218] = 1.462694e-024;
m[-1218] = 1.462694e-024;
m[2101] = 0.000000e+000;
m[-2101] = 0.000000e+000;
m[2103] = 0.000000e+000;
m[-2103] = 0.000000e+000;
m[2116] = 4.388081e-024;
m[-2116] = 4.388081e-024;
m[2118] = 2.194041e-024;
m[-2118] = 2.194041e-024;
m[2122] = 4.388081e-024;
m[-2122] = 4.388081e-024;
m[2124] = 5.485102e-024;
m[-2124] = 5.485102e-024;
m[2126] = 1.880606e-024;
m[-2126] = 1.880606e-024;
m[2128] = 1.462694e-024;
m[-2128] = 1.462694e-024;
m[2203] = 0.000000e+000;
m[-2203] = 0.000000e+000;
m[2216] = 4.388081e-024;
m[-2216] = 4.388081e-024;
m[2218] = 2.194041e-024;
m[-2218] = 2.194041e-024;
m[2222] = 4.388081e-024;
m[-2222] = 4.388081e-024;
m[2226] = 1.880606e-024;
m[-2226] = 1.880606e-024;
m[2228] = 2.194041e-024;
m[-2228] = 2.194041e-024;
m[3101] = 0.000000e+000;
m[-3101] = 0.000000e+000;
m[3103] = 0.000000e+000;
m[-3103] = 0.000000e+000;
m[3114] = 1.670589e-023;
m[-3114] = 1.670589e-023;
m[3116] = 5.485102e-024;
m[-3116] = 5.485102e-024;
m[3118] = 3.656734e-024;
m[-3118] = 3.656734e-024;
m[3124] = 4.219309e-023;
m[-3124] = 4.219309e-023;
m[3126] = 8.227653e-024;
m[-3126] = 8.227653e-024;
m[3128] = 3.291061e-024;
m[-3128] = 3.291061e-024;
m[3201] = 0.000000e+000;
m[-3201] = 0.000000e+000;
m[3203] = 0.000000e+000;
m[-3203] = 0.000000e+000;
m[3214] = 1.828367e-023;
m[-3214] = 1.828367e-023;
m[3216] = 5.485102e-024;
m[-3216] = 5.485102e-024;
m[3218] = 3.656734e-024;
m[-3218] = 3.656734e-024;
m[3224] = 1.838582e-023;
m[-3224] = 1.838582e-023;
m[3226] = 5.485102e-024;
m[-3226] = 5.485102e-024;
m[3228] = 3.656734e-024;
m[-3228] = 3.656734e-024;
m[3303] = 0.000000e+000;
m[-3303] = 0.000000e+000;
m[3314] = 6.648608e-023;
m[-3314] = 6.648608e-023;
m[3324] = 7.233101e-023;
m[-3324] = 7.233101e-023;
m[4101] = 0.000000e+000;
m[-4101] = 0.000000e+000;
m[4103] = 0.000000e+000;
m[-4103] = 0.000000e+000;
m[4114] = 0.000000e+000;
m[-4114] = 0.000000e+000;
m[4201] = 0.000000e+000;
m[-4201] = 0.000000e+000;
m[4203] = 0.000000e+000;
m[-4203] = 0.000000e+000;
m[4214] = 3.291061e-022;
m[-4214] = 3.291061e-022;
m[4224] = 0.000000e+000;
m[-4224] = 0.000000e+000;
m[4301] = 0.000000e+000;
m[-4301] = 0.000000e+000;
m[4303] = 0.000000e+000;
m[-4303] = 0.000000e+000;
m[4314] = 0.000000e+000;
m[-4314] = 0.000000e+000;
m[4324] = 0.000000e+000;
m[-4324] = 0.000000e+000;
m[4334] = 0.000000e+000;
m[-4334] = 0.000000e+000;
m[4403] = 0.000000e+000;
m[-4403] = 0.000000e+000;
m[4412] = 3.335641e-013;
m[-4412] = 3.335641e-013;
m[4414] = 3.335641e-013;
m[-4414] = 3.335641e-013;
m[4422] = 3.335641e-013;
m[-4422] = 3.335641e-013;
m[4424] = 3.335641e-013;
m[-4424] = 3.335641e-013;
m[4432] = 3.335641e-013;
m[-4432] = 3.335641e-013;
m[4434] = 3.335641e-013;
m[-4434] = 3.335641e-013;
m[4444] = 3.335641e-013;
m[-4444] = 3.335641e-013;
m[5101] = 0.000000e+000;
m[-5101] = 0.000000e+000;
m[5103] = 0.000000e+000;
m[-5103] = 0.000000e+000;
m[5114] = 0.000000e+000;
m[-5114] = 0.000000e+000;
m[5142] = 1.290893e-012;
m[-5142] = 1.290893e-012;
m[5201] = 0.000000e+000;
m[-5201] = 0.000000e+000;
m[5203] = 0.000000e+000;
m[-5203] = 0.000000e+000;
m[5214] = 0.000000e+000;
m[-5214] = 0.000000e+000;
m[5224] = 0.000000e+000;
m[-5224] = 0.000000e+000;
m[5242] = 1.290893e-012;
m[-5242] = 1.290893e-012;
m[5301] = 0.000000e+000;
m[-5301] = 0.000000e+000;
m[5303] = 0.000000e+000;
m[-5303] = 0.000000e+000;
m[5314] = 0.000000e+000;
m[-5314] = 0.000000e+000;
m[5324] = 0.000000e+000;
m[-5324] = 0.000000e+000;
m[5334] = 0.000000e+000;
m[-5334] = 0.000000e+000;
m[5342] = 1.290893e-012;
m[-5342] = 1.290893e-012;
m[5401] = 0.000000e+000;
m[-5401] = 0.000000e+000;
m[5403] = 0.000000e+000;
m[-5403] = 0.000000e+000;
m[5412] = 1.290893e-012;
m[-5412] = 1.290893e-012;
m[5414] = 1.290893e-012;
m[-5414] = 1.290893e-012;
m[5422] = 1.290893e-012;
m[-5422] = 1.290893e-012;
m[5424] = 1.290893e-012;
m[-5424] = 1.290893e-012;
m[5432] = 1.290893e-012;
m[-5432] = 1.290893e-012;
m[5434] = 1.290893e-012;
m[-5434] = 1.290893e-012;
m[5442] = 1.290893e-012;
m[-5442] = 1.290893e-012;
m[5444] = 1.290893e-012;
m[-5444] = 1.290893e-012;
m[5503] = 0.000000e+000;
m[-5503] = 0.000000e+000;
m[5512] = 1.290893e-012;
m[-5512] = 1.290893e-012;
m[5514] = 1.290893e-012;
m[-5514] = 1.290893e-012;
m[5522] = 1.290893e-012;
m[-5522] = 1.290893e-012;
m[5524] = 1.290893e-012;
m[-5524] = 1.290893e-012;
m[5532] = 1.290893e-012;
m[-5532] = 1.290893e-012;
m[5534] = 1.290893e-012;
m[-5534] = 1.290893e-012;
m[5542] = 1.290893e-012;
m[-5542] = 1.290893e-012;
m[5544] = 1.290893e-012;
m[-5544] = 1.290893e-012;
m[5554] = 1.290893e-012;
m[-5554] = 1.290893e-012;
m[10022] = 0.000000e+000;
m[10111] = 2.483820e-024;
m[10113] = 4.635297e-024;
m[10115] = 2.541360e-024;
m[10211] = 2.483820e-024;
m[-10211] = 2.483820e-024;
m[10213] = 4.635297e-024;
m[-10213] = 4.635297e-024;
m[10215] = 2.541360e-024;
m[-10215] = 2.541360e-024;
m[9010221] = 1.316424e-023;
m[10223] = 1.828367e-024;
m[10225] = 0.000000e+000;
m[10315] = 3.538775e-024;
m[-10315] = 3.538775e-024;
m[10325] = 3.538775e-024;
m[-10325] = 3.538775e-024;
m[10331] = 5.265698e-024;
m[10333] = 0.000000e+000;
m[10335] = 0.000000e+000;
m[10443] = 0.000000e+000;
m[10541] = 0.000000e+000;
m[-10541] = 0.000000e+000;
m[10543] = 0.000000e+000;
m[-10543] = 0.000000e+000;
m[10551] = 1.000000e+016;
m[10553] = 0.000000e+000;
m[10555] = 0.000000e+000;
m[11112] = 0.000000e+000;
m[-11112] = 0.000000e+000;
m[11114] = 2.194041e-024;
m[-11114] = 2.194041e-024;
m[11116] = 1.880606e-024;
m[-11116] = 1.880606e-024;
m[11212] = 1.880606e-024;
m[-11212] = 1.880606e-024;
m[11216] = 0.000000e+000;
m[-11216] = 0.000000e+000;
m[12112] = 1.880606e-024;
m[-12112] = 1.880606e-024;
m[12114] = 2.194041e-024;
m[-12114] = 2.194041e-024;
m[12116] = 5.063171e-024;
m[-12116] = 5.063171e-024;
m[12118] = 0.000000e+000;
m[-12118] = 0.000000e+000;
m[12122] = 0.000000e+000;
m[-12122] = 0.000000e+000;
m[12126] = 1.880606e-024;
m[-12126] = 1.880606e-024;
m[12212] = 1.880606e-024;
m[-12212] = 1.880606e-024;
m[12214] = 2.194041e-024;
m[-12214] = 2.194041e-024;
m[12216] = 5.063171e-024;
m[-12216] = 5.063171e-024;
m[12218] = 0.000000e+000;
m[-12218] = 0.000000e+000;
m[12222] = 0.000000e+000;
m[-12222] = 0.000000e+000;
m[12224] = 2.194041e-024;
m[-12224] = 2.194041e-024;
m[12226] = 1.880606e-024;
m[-12226] = 1.880606e-024;
m[13112] = 6.582122e-024;
m[-13112] = 6.582122e-024;
m[13114] = 1.097020e-023;
m[-13114] = 1.097020e-023;
m[13116] = 5.485102e-024;
m[-13116] = 5.485102e-024;
m[13122] = 1.316424e-023;
m[-13122] = 1.316424e-023;
m[13124] = 1.097020e-023;
m[-13124] = 1.097020e-023;
m[13126] = 6.928549e-024;
m[-13126] = 6.928549e-024;
m[13212] = 6.582122e-024;
m[-13212] = 6.582122e-024;
m[13214] = 1.097020e-023;
m[-13214] = 1.097020e-023;
m[13216] = 5.485102e-024;
m[-13216] = 5.485102e-024;
m[13222] = 6.582122e-024;
m[-13222] = 6.582122e-024;
m[13224] = 1.097020e-023;
m[-13224] = 1.097020e-023;
m[13226] = 5.485102e-024;
m[-13226] = 5.485102e-024;
m[13314] = 2.742551e-023;
m[-13314] = 2.742551e-023;
m[13316] = 0.000000e+000;
m[-13316] = 0.000000e+000;
m[13324] = 2.742551e-023;
m[-13324] = 2.742551e-023;
m[13326] = 0.000000e+000;
m[-13326] = 0.000000e+000;
m[14122] = 1.828367e-022;
m[-14122] = 1.828367e-022;
m[14124] = 0.000000e+000;
m[-14124] = 0.000000e+000;
m[10221] = 2.194040e-024;
m[20223] = 2.742551e-023;
m[20315] = 2.384827e-024;
m[-20315] = 2.384827e-024;
m[20325] = 2.384827e-024;
m[-20325] = 2.384827e-024;
m[20333] = 1.185968e-023;
m[20543] = 0.000000e+000;
m[-20543] = 0.000000e+000;
m[20553] = 1.000000e+016;
m[20555] = 0.000000e+000;
m[21112] = 2.632849e-024;
m[-21112] = 2.632849e-024;
m[21114] = 3.291061e-024;
m[-21114] = 3.291061e-024;
m[21212] = 2.632849e-024;
m[-21212] = 2.632849e-024;
m[21214] = 6.582122e-024;
m[-21214] = 6.582122e-024;
m[22112] = 4.388081e-024;
m[-22112] = 4.388081e-024;
m[22114] = 3.291061e-024;
m[-22114] = 3.291061e-024;
m[22122] = 2.632849e-024;
m[-22122] = 2.632849e-024;
m[22124] = 6.582122e-024;
m[-22124] = 6.582122e-024;
m[22212] = 4.388081e-024;
m[-22212] = 4.388081e-024;
m[22214] = 3.291061e-024;
m[-22214] = 3.291061e-024;
m[22222] = 2.632849e-024;
m[-22222] = 2.632849e-024;
m[22224] = 3.291061e-024;
m[-22224] = 3.291061e-024;
m[23112] = 7.313469e-024;
m[-23112] = 7.313469e-024;
m[23114] = 2.991874e-024;
m[-23114] = 2.991874e-024;
m[23122] = 4.388081e-024;
m[-23122] = 4.388081e-024;
m[23124] = 6.582122e-024;
m[-23124] = 6.582122e-024;
m[23126] = 3.291061e-024;
m[-23126] = 3.291061e-024;
m[23212] = 7.313469e-024;
m[-23212] = 7.313469e-024;
m[23214] = 2.991874e-024;
m[-23214] = 2.991874e-024;
m[23222] = 7.313469e-024;
m[-23222] = 7.313469e-024;
m[23224] = 2.991874e-024;
m[-23224] = 2.991874e-024;
m[23314] = 0.000000e+000;
m[-23314] = 0.000000e+000;
m[23324] = 0.000000e+000;
m[-23324] = 0.000000e+000;
m[30113] = 2.742551e-024;
m[30213] = 2.742551e-024;
m[-30213] = 2.742551e-024;
m[30223] = 2.991874e-024;
m[30313] = 2.056913e-024;
m[-30313] = 2.056913e-024;
m[30323] = 2.056913e-024;
m[-30323] = 2.056913e-024;
m[30343] = 0.000000e+000;
m[-30343] = 0.000000e+000;
m[30353] = 0.000000e+000;
m[-30353] = 0.000000e+000;
m[30363] = 0.000000e+000;
m[-30363] = 0.000000e+000;
m[30411] = 0.000000e+000;
m[-30411] = 0.000000e+000;
m[30413] = 0.000000e+000;
m[-30413] = 0.000000e+000;
m[30421] = 0.000000e+000;
m[-30421] = 0.000000e+000;
m[30423] = 0.000000e+000;
m[-30423] = 0.000000e+000;
m[30443] = 2.789035e-023;
m[30553] = 0.000000e+000;
m[31114] = 1.880606e-024;
m[-31114] = 1.880606e-024;
m[31214] = 4.388081e-024;
m[-31214] = 4.388081e-024;
m[32112] = 4.388081e-024;
m[-32112] = 4.388081e-024;
m[32114] = 1.880606e-024;
m[-32114] = 1.880606e-024;
m[32124] = 4.388081e-024;
m[-32124] = 4.388081e-024;
m[32212] = 4.388081e-024;
m[-32212] = 4.388081e-024;
m[32214] = 1.880606e-024;
m[-32214] = 1.880606e-024;
m[32224] = 1.880606e-024;
m[-32224] = 1.880606e-024;
m[33122] = 1.880606e-023;
m[-33122] = 1.880606e-023;
m[33314] = 0.000000e+000;
m[-33314] = 0.000000e+000;
m[33324] = 0.000000e+000;
m[-33324] = 0.000000e+000;
m[41214] = 0.000000e+000;
m[-41214] = 0.000000e+000;
m[42112] = 6.582122e-024;
m[-42112] = 6.582122e-024;
m[42124] = 0.000000e+000;
m[-42124] = 0.000000e+000;
m[42212] = 6.582122e-024;
m[-42212] = 6.582122e-024;
m[43122] = 2.194041e-024;
m[-43122] = 2.194041e-024;
m[52114] = 0.000000e+000;
m[-52114] = 0.000000e+000;
m[52214] = 0.000000e+000;
m[-52214] = 0.000000e+000;
m[53122] = 4.388081e-024;
m[-53122] = 4.388081e-024;
m[100111] = 1.645531e-024;
m[100113] = 2.123265e-024;
m[100211] = 1.645531e-024;
m[-100211] = 1.645531e-024;
m[100213] = 2.123265e-024;
m[-100213] = 2.123265e-024;
m[100221] = 1.196749e-023;
m[100223] = 3.871836e-024;
m[100225] = 0.000000e+000;
m[100311] = 0.000000e+000;
m[-100311] = 0.000000e+000;
m[100313] = 2.837122e-024;
m[-100313] = 2.837122e-024;
m[100315] = 0.000000e+000;
m[-100315] = 0.000000e+000;
m[100321] = 0.000000e+000;
m[-100321] = 0.000000e+000;
m[100323] = 2.837122e-024;
m[-100323] = 2.837122e-024;
m[100325] = 0.000000e+000;
m[-100325] = 0.000000e+000;
m[100331] = 0.000000e+000;
m[100333] = 4.388081e-024;
m[100335] = 3.291061e-024;
m[100441] = 0.000000e+000;
m[100551] = 0.000000e+000;
m[100553] = 1.495937e-020;
m[100555] = 1.000000e+016;
m[100557] = 0.000000e+000;
m[110551] = 1.000000e+016;
m[110553] = 0.000000e+000;
m[110555] = 0.000000e+000;
m[120553] = 1.000000e+016;
m[120555] = 0.000000e+000;
m[130553] = 0.000000e+000;
m[200111] = 3.134344e-024;
m[200211] = 3.134344e-024;
m[-200211] = 3.134344e-024;
m[200551] = 0.000000e+000;
m[200553] = 2.502708e-020;
m[200555] = 0.000000e+000;
m[210551] = 0.000000e+000;
m[210553] = 0.000000e+000;
m[220553] = 0.000000e+000;
m[300553] = 4.701516e-023;
m[9000221] = 0.000000e+000;
m[9000443] = 1.265793e-023;
m[9000553] = 5.983747e-024;
m[9010443] = 8.438618e-024;
m[9010553] = 8.331800e-024;
m[9020221] = 6.038644e-024;
m[9020443] = 1.530726e-023;
m[9060225] = 4.388081e-024;
m[9070225] = 2.056913e-024;
m[1000001] = 0.000000e+000;
m[-1000001] = 0.000000e+000;
m[1000002] = 0.000000e+000;
m[-1000002] = 0.000000e+000;
m[1000003] = 0.000000e+000;
m[-1000003] = 0.000000e+000;
m[1000004] = 0.000000e+000;
m[-1000004] = 0.000000e+000;
m[1000005] = 0.000000e+000;
m[-1000005] = 0.000000e+000;
m[1000006] = 0.000000e+000;
m[-1000006] = 0.000000e+000;
m[1000011] = 0.000000e+000;
m[-1000011] = 0.000000e+000;
m[1000012] = 0.000000e+000;
m[-1000012] = 0.000000e+000;
m[1000013] = 0.000000e+000;
m[-1000013] = 0.000000e+000;
m[1000014] = 0.000000e+000;
m[-1000014] = 0.000000e+000;
m[1000015] = 0.000000e+000;
m[-1000015] = 0.000000e+000;
m[1000016] = 0.000000e+000;
m[-1000016] = 0.000000e+000;
m[1000021] = 0.000000e+000;
m[1000022] = 0.000000e+000;
m[1000023] = 0.000000e+000;
m[1000024] = 0.000000e+000;
m[-1000024] = 0.000000e+000;
m[1000025] = 0.000000e+000;
m[1000035] = 0.000000e+000;
m[1000037] = 0.000000e+000;
m[-1000037] = 0.000000e+000;
m[1000039] = 0.000000e+000;
m[2000001] = 0.000000e+000;
m[-2000001] = 0.000000e+000;
m[2000002] = 0.000000e+000;
m[-2000002] = 0.000000e+000;
m[2000003] = 0.000000e+000;
m[-2000003] = 0.000000e+000;
m[2000004] = 0.000000e+000;
m[-2000004] = 0.000000e+000;
m[2000005] = 0.000000e+000;
m[-2000005] = 0.000000e+000;
m[2000006] = 0.000000e+000;
m[-2000006] = 0.000000e+000;
m[2000011] = 0.000000e+000;
m[-2000011] = 0.000000e+000;
m[2000012] = 0.000000e+000;
m[-2000012] = 0.000000e+000;
m[2000013] = 0.000000e+000;
m[-2000013] = 0.000000e+000;
m[2000014] = 0.000000e+000;
m[-2000014] = 0.000000e+000;
m[2000015] = 0.000000e+000;
m[-2000015] = 0.000000e+000;
m[2000016] = 0.000000e+000;
m[-2000016] = 0.000000e+000;
m[3000111] = 0.000000e+000;
m[3000113] = 0.000000e+000;
m[3000211] = 0.000000e+000;
m[-3000211] = 0.000000e+000;
m[3000213] = 0.000000e+000;
m[-3000213] = 0.000000e+000;
m[3000221] = 0.000000e+000;
m[3000223] = 0.000000e+000;
m[3000331] = 0.000000e+000;
m[3100021] = 0.000000e+000;
m[3100111] = 0.000000e+000;
m[3100113] = 0.000000e+000;
m[3200111] = 0.000000e+000;
m[3200113] = 0.000000e+000;
m[3300113] = 0.000000e+000;
m[3400113] = 0.000000e+000;
m[4000001] = 0.000000e+000;
m[-4000001] = 0.000000e+000;
m[4000002] = 0.000000e+000;
m[-4000002] = 0.000000e+000;
m[4000011] = 0.000000e+000;
m[-4000011] = 0.000000e+000;
m[4000012] = 0.000000e+000;
m[-4000012] = 0.000000e+000;
m[5000039] = 0.000000e+000;
m[9900012] = 0.000000e+000;
m[9900014] = 0.000000e+000;
m[9900016] = 0.000000e+000;
m[9900023] = 0.000000e+000;
m[9900024] = 0.000000e+000;
m[-9900024] = 0.000000e+000;
m[9900041] = 0.000000e+000;
m[-9900041] = 0.000000e+000;
m[9900042] = 0.000000e+000;
m[-9900042] = 0.000000e+000;
m[1027013000] = 0.000000e+000;
m[1012006000] = 0.000000e+000;
m[1063029000] = 0.000000e+000;
m[1014007000] = 0.000000e+000;
m[1016008000] = 0.000000e+000;
m[1028014000] = 0.000000e+000;
m[1065029000] = 0.000000e+000;
m[1009004000] = 0.000000e+000;
m[1019009000] = 0.000000e+000;
m[1056026000] = 0.000000e+000;
m[1207082000] = 0.000000e+000;
m[1208082000] = 0.000000e+000;
m[1029014000] = 0.000000e+000;
m[1206082000] = 0.000000e+000;
m[1054026000] = 0.000000e+000;
m[1018008000] = 0.000000e+000;
m[1030014000] = 0.000000e+000;
m[1057026000] = 0.000000e+000;
m[1204082000] = 0.000000e+000;
m[-99000000] = 0.000000e+000;
m[1028013000] = 0.000000e+000;
m[1040018000] = 0.000000e+000;
m[1011005000] = 0.000000e+000;
m[1012005000] = 0.000000e+000;
m[1013006000] = 0.000000e+000;
m[1014006000] = 0.000000e+000;
m[1052024000] = 0.000000e+000;
m[1024012000] = 0.000000e+000;
m[1026012000] = 0.000000e+000;
m[1027012000] = 0.000000e+000;
m[1015007000] = 0.000000e+000;
m[1022010000] = 0.000000e+000;
m[1058028000] = 0.000000e+000;
m[1060028000] = 0.000000e+000;
m[1062028000] = 0.000000e+000;
m[1064028000] = 0.000000e+000;
m[1007003000] = 0.000000e+000;
m[1025012000] = 0.000000e+000;
m[1053024000] = 0.000000e+000;
m[1055025000] = 0.000000e+000;
m[1008004000] = 0.000000e+000;
m[1010004000] = 0.000000e+000;
m[1010005000] = 0.000000e+000;
m[1016007000] = 0.000000e+000;
m[1017008000] = 0.000000e+000;
m[1019008000] = 0.000000e+000;
m[1023010000] = 0.000000e+000;
m[1024011000] = 0.000000e+000;
m[1031015000] = 0.000000e+000;
m[1039017000] = 0.000000e+000;
m[1040017000] = 0.000000e+000;
m[1036018000] = 0.000000e+000;
m[1050024000] = 0.000000e+000;
m[1054024000] = 0.000000e+000;
m[1059026000] = 0.000000e+000;
m[1061028000] = 0.000000e+000;
m[1063028000] = 0.000000e+000;
m[1092042000] = 0.000000e+000;
m[1095042000] = 0.000000e+000;
m[1096042000] = 0.000000e+000;
m[1097042000] = 0.000000e+000;
m[1098042000] = 0.000000e+000;
m[1100042000] = 0.000000e+000;
m[1108046000] = 0.000000e+000;
// Added by hand:
m[9902210] = 0.000000e+000; //diffractive p-state -> assume no lifetime
return true;
}
private:
/// @name Histograms
//@{
Histo1DPtr _h_mult_total; // full kinematic range
Histo1DPtr _h_mult_eta[5]; // in eta bins
Histo1DPtr _h_mult_pt[5]; // in pT bins
Histo1DPtr _h_dndeta; // density dn/deta
Histo1DPtr _h_dndpt; // density dn/dpT
//@}
/// @name Private variables
double _p_min;
double _pt_min;
double _eta_min;
double _eta_max;
double _maxlft;
/// Count selected events
- double _sumW;
+ CounterPtr _sumW;
map<int, double> _partLftMap; // Map <PDGID, PDGLIFETIME>
};
// The hook for the plugin system
DECLARE_RIVET_PLUGIN(LHCB_2014_I1281685);
}
diff --git a/analyses/pluginLHCb/LHCB_2015_I1333223.cc b/analyses/pluginLHCb/LHCB_2015_I1333223.cc
--- a/analyses/pluginLHCb/LHCB_2015_I1333223.cc
+++ b/analyses/pluginLHCb/LHCB_2015_I1333223.cc
@@ -1,113 +1,112 @@
// -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Tools/Logging.hh"
#include "Rivet/Projections/ChargedFinalState.hh"
#include "Rivet/Math/Units.hh"
#include <vector>
using namespace std;
namespace Rivet {
class LHCB_2015_I1333223 : public Analysis {
public:
/// @name Constructors etc.
//@{
/// Constructor
LHCB_2015_I1333223()
: Analysis("LHCB_2015_I1333223")
{ }
//@}
public:
/// @name Analysis methods
//@{
/// Book histograms and initialise projections before the run
void init() {
// Charged particles
declare(ChargedFinalState(Cuts::eta> 2.0 && Cuts::eta <4.5 && Cuts::pT >0.2*GeV), "CFS");
// Reproducing only measurement for prompt charged particles
book(_hInelasticXs ,1, 1, 1);
}
/// Perform the per-event analysis
void analyze(const Event& event) {
- const double weight = 1.0;
const ChargedFinalState &cfs = apply<ChargedFinalState> (event, "CFS");
// eliminate non-inelastic events and empty events in LHCb
if (cfs.particles().size() == 0) vetoEvent;
// See if this event has at least one prompt particle
foreach (const Particle &myp, cfs.particles()){
double dPV = getPVDCA(myp);
// if IP > 200 microns the particle is not considered prompt
if ((dPV < 0.) || (dPV > 0.2 * millimeter)) {
MSG_DEBUG(" Vetoing " << myp.pdgId() << " at " << dPV);
continue;
}
// histo gets filled only for inelastic events (at least one prompt charged particle)
- _hInelasticXs->fill(sqrtS(), weight);
+ _hInelasticXs->fill(sqrtS());
break;
} //end loop on particles
}
/// Normalise histograms etc., after the run
void finalize() {
scale(_hInelasticXs, crossSection()/sumOfWeights()/millibarn);
}
//@}
private:
/*
* Compute Distance of Closest Approach in z range for one particle.
* Assuming length unit is mm.
* Returns -1. if unable to compute the DCA to PV.
*/
double getPVDCA(const Particle& p) {
const HepMC::GenVertex* vtx = p.genParticle()->production_vertex();
if ( 0 == vtx ) return -1.;
// Unit vector of particle's MOMENTUM three vector
const Vector3 u = p.momentum().p3().unit();
// The interaction point is always at (0, 0,0,0) hence the
// vector pointing from the PV to the particle production vertex is:
Vector3 d(vtx->position().x(), vtx->position().y(), vtx->position().z());
// Subtract projection of d onto u from d
double proj = d.dot(u);
d -= (u * proj);
// d should be orthogonal to u and it's length give the distance of
// closest approach
return d.mod();
}
/// @name Histograms
//@{
Histo1DPtr _hInelasticXs;
//@}
//
};
// The hook for the plugin system
DECLARE_RIVET_PLUGIN(LHCB_2015_I1333223);
}

File Metadata

Mime Type
text/x-diff
Expires
Tue, Nov 19, 7:12 PM (1 d, 12 h)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
3805778
Default Alt Text
(127 KB)

Event Timeline