diff --git a/analyses/pluginBES/BESIII_2012_I1113599.cc b/analyses/pluginBES/BESIII_2012_I1113599.cc --- a/analyses/pluginBES/BESIII_2012_I1113599.cc +++ b/analyses/pluginBES/BESIII_2012_I1113599.cc @@ -1,126 +1,130 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/Beam.hh" #include "Rivet/Projections/FinalState.hh" namespace Rivet { /// @brief J/psi to p pbar n nbar class BESIII_2012_I1113599 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2012_I1113599); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(Beam(), "Beams"); declare(UnstableParticles(), "UFS"); declare(FinalState(), "FS"); // Book histograms - _h_proton = bookHisto1D("ctheta_p",20,-1.,1.); - _h_neutron = bookHisto1D("ctheta_n",20,-1.,1.); + book(_h_proton , "ctheta_p",20,-1.,1.); + book(_h_neutron, "ctheta_n",20,-1.,1.); } /// Perform the per-event analysis void analyze(const Event& event) { // get the axis, direction of incoming electron const ParticlePair& beams = apply(event, "Beams").beams(); Vector3 axis; - if(beams.first.pdgId()>0) + if(beams.first.pid()>0) axis = beams.first .momentum().p3().unit(); else axis = beams.second.momentum().p3().unit(); // types of final state particles const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); Particle outgoing; - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; - if(p.pdgId()==2212 || p.pdgId()==2112) + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; + if(p.pid()==2212 || p.pid()==2112) outgoing = p; ++ntotal; } if(ntotal==2) { if(nCount[2212]==1 && nCount[-2212]==1) { - _h_proton->fill(outgoing.momentum().p3().unit().dot(axis),event.weight()); + _h_proton->fill(outgoing.momentum().p3().unit().dot(axis)); } else if(nCount[2112]==1 && nCount[-2112]==1) { - _h_neutron->fill(outgoing.momentum().p3().unit().dot(axis),event.weight()); + _h_neutron->fill(outgoing.momentum().p3().unit().dot(axis)); } } } pair > calcAlpha(Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,make_pair(0.,0.)); double sum1(0.),sum2(0.),sum3(0.),sum4(0.),sum5(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double a = 1.5*(bin.xMax() - bin.xMin()); double b = 0.5*(pow(bin.xMax(),3) - pow(bin.xMin(),3)); double Ei = bin.areaErr(); sum1 += a*Oi/sqr(Ei); sum2 += b*Oi/sqr(Ei); sum3 += sqr(a)/sqr(Ei); sum4 += sqr(b)/sqr(Ei); sum5 += a*b/sqr(Ei); } // calculate alpha double alpha = (-3*sum1 + 9*sum2 + sum3 - 3*sum5)/(sum1 - 3*sum2 + 3*sum4 - sum5); // and error double cc = -pow((sum3 + 9*sum4 - 6*sum5),3); double bb = -2*sqr(sum3 + 9*sum4 - 6*sum5)*(sum1 - 3*sum2 + 3*sum4 - sum5); double aa = sqr(sum1 - 3*sum2 + 3*sum4 - sum5)*(-sum3 - 9*sum4 + sqr(sum1 - 3*sum2 + 3*sum4 - sum5) + 6*sum5); double dis = sqr(bb)-4.*aa*cc; if(dis>0.) { dis = sqrt(dis); return make_pair(alpha,make_pair(0.5*(-bb+dis)/aa,-0.5*(-bb-dis)/aa)); } else { return make_pair(alpha,make_pair(0.,0.)); } } /// Normalise histograms etc., after the run void finalize() { // proton normalize(_h_proton ); pair > alpha = calcAlpha(_h_proton); - Scatter2DPtr _h_alpha_proton = bookScatter2D(1,1,1); - _h_alpha_proton->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second.first,alpha.second.second) ); + Scatter2DPtr _h_alpha_proton; + book(_h_alpha_proton,1,1,1); + _h_alpha_proton->addPoint(0.5, alpha.first, make_pair(0.5,0.5), + make_pair(alpha.second.first,alpha.second.second) ); // neutron normalize(_h_neutron); alpha = calcAlpha(_h_neutron); - Scatter2DPtr _h_alpha_neutron = bookScatter2D(1,1,2); - _h_alpha_neutron->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second.first,alpha.second.second) ); + Scatter2DPtr _h_alpha_neutron; + book(_h_alpha_neutron, 1,1,2); + _h_alpha_neutron->addPoint(0.5, alpha.first, make_pair(0.5,0.5), + make_pair(alpha.second.first,alpha.second.second) ); } //@} /// @name Histograms //@{ Histo1DPtr _h_proton,_h_neutron; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2012_I1113599); } diff --git a/analyses/pluginBES/BESIII_2014_I1323621.cc b/analyses/pluginBES/BESIII_2014_I1323621.cc --- a/analyses/pluginBES/BESIII_2014_I1323621.cc +++ b/analyses/pluginBES/BESIII_2014_I1323621.cc @@ -1,114 +1,115 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2014_I1323621 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2014_I1323621); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); - _nChi0 = bookCounter("TMP/chi0"); + book(_nChi0, "TMP/chi0"); } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - --nRes[child.pdgId()]; + --nRes[child.pid()]; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } const FinalState& ufs = apply(event, "UFS"); bool matched = false; - foreach (const Particle& p, ufs.particles(Cuts::pid==10441)) { + for (const Particle& p : ufs.particles(Cuts::pid==10441)) { if(p.children().empty()) continue; map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); - foreach (const Particle & p2, ufs.particles(Cuts::pid==223)) { + for (const Particle& p2 : ufs.particles(Cuts::pid==223)) { map nRes2 = nRes; int ncount2 = ncount; findChildren(p2,nRes2,ncount2); if(ncount2!=0) continue; matched=true; for(auto const & val : nRes2) { if(val.second!=0) { matched = false; break; } } if(matched) { - _nChi0->fill(event.weight()); + _nChi0->fill(); break; } } if(matched) break; } } /// Normalise histograms etc., after the run void finalize() { double fact = crossSection()/ sumOfWeights() /picobarn; double sigma = _nChi0->val()*fact; double error = _nChi0->err()*fact; Scatter2D temphisto(refData(1, 1, 1)); - Scatter2DPtr mult = bookScatter2D(1, 1, 1); + Scatter2DPtr mult; + book(mult, 1, 1, 1); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/GeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } //@} /// @name Histograms //@{ CounterPtr _nChi0; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2014_I1323621); } diff --git a/analyses/pluginBES/BESIII_2015_I1329785.cc b/analyses/pluginBES/BESIII_2015_I1329785.cc --- a/analyses/pluginBES/BESIII_2015_I1329785.cc +++ b/analyses/pluginBES/BESIII_2015_I1329785.cc @@ -1,132 +1,133 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2015_I1329785 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2015_I1329785); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); - _nChi0 = bookCounter("TMP/chi0"); - _nChi1 = bookCounter("TMP/chi1"); - _nChi2 = bookCounter("TMP/chi2"); + book(_nChi0, "TMP/chi0"); + book(_nChi1, "TMP/chi1"); + book(_nChi2, "TMP/chi2"); } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - --nRes[child.pdgId()]; + --nRes[child.pid()]; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } const FinalState& ufs = apply(event, "UFS"); - foreach (const Particle& p, ufs.particles(Cuts::pid==10441 or Cuts::pid==20443 or Cuts::pid==445)) { + for (const Particle& p : ufs.particles(Cuts::pid==10441 or Cuts::pid==20443 or Cuts::pid==445)) { if(p.children().empty()) continue; map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); // chi gamma if(ncount!=1) continue; bool matched = true; for(auto const & val : nRes) { if(val.first==22) { if(val.second !=1) { matched = false; break; } } else if(val.second!=0) { matched = false; break; } } if(matched) { - if(p.pdgId()==10441) - _nChi0->fill(event.weight()); - else if(p.pdgId()==20443) - _nChi1->fill(event.weight()); - else if(p.pdgId()==445) - _nChi2->fill(event.weight()); + if(p.pid()==10441) + _nChi0->fill(); + else if(p.pid()==20443) + _nChi1->fill(); + else if(p.pid()==445) + _nChi2->fill(); break; } } } /// Normalise histograms etc., after the run void finalize() { double fact = crossSection()/ sumOfWeights() /picobarn; for(unsigned int ix=1;ix<4;++ix) { double sigma(0.),error(0.); if(ix==1) { sigma = _nChi0->val()*fact; error = _nChi0->err()*fact; } else if(ix==2) { sigma = _nChi1->val()*fact; error = _nChi1->err()*fact; } else if(ix==3) { sigma = _nChi2->val()*fact; error = _nChi2->err()*fact; } Scatter2D temphisto(refData(ix, 1, 8)); - Scatter2DPtr mult = bookScatter2D(ix, 1, 8); + Scatter2DPtr mult; + book(mult,ix, 1, 8); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/GeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } } //@} /// @name Histograms //@{ CounterPtr _nChi0,_nChi1,_nChi2; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2015_I1329785); } diff --git a/analyses/pluginBES/BESIII_2015_I1406939.cc b/analyses/pluginBES/BESIII_2015_I1406939.cc --- a/analyses/pluginBES/BESIII_2015_I1406939.cc +++ b/analyses/pluginBES/BESIII_2015_I1406939.cc @@ -1,129 +1,130 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2015_I1406939 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2015_I1406939); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); - _nChi1 = bookCounter("TMP/chi1"); - _nChi2 = bookCounter("TMP/chi2"); + book(_nChi1, "TMP/chi1"); + book(_nChi2, "TMP/chi2"); } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - --nRes[child.pdgId()]; + --nRes[child.pid()]; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } const FinalState& ufs = apply(event, "UFS"); bool matched = false; - foreach (const Particle& p, ufs.particles(Cuts::pid==20443 or Cuts::pid==445)) { + for (const Particle& p : ufs.particles(Cuts::pid==20443 or Cuts::pid==445)) { if(p.children().empty()) continue; map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); - foreach (const Particle & p2, ufs.particles(Cuts::pid==223)) { + for (const Particle& p2 : ufs.particles(Cuts::pid==223)) { map nRes2 = nRes; int ncount2 = ncount; findChildren(p2,nRes2,ncount2); if(ncount2!=0) continue; matched=true; for(auto const & val : nRes2) { if(val.second!=0) { matched = false; break; } } if(matched) { - if(p.pdgId()==20443) - _nChi1->fill(event.weight()); - else if(p.pdgId()==445) - _nChi2->fill(event.weight()); + if(p.pid()==20443) + _nChi1->fill(); + else if(p.pid()==445) + _nChi2->fill(); break; } } if(matched) break; } } /// Normalise histograms etc., after the run void finalize() { double fact = crossSection()/ sumOfWeights() /picobarn; for(unsigned int ix=2;ix<4;++ix) { double sigma(0.),error(0.); if(ix==2) { sigma = _nChi1->val()*fact; error = _nChi1->err()*fact; } else if(ix==3) { sigma = _nChi2->val()*fact; error = _nChi2->err()*fact; } Scatter2D temphisto(refData(ix, 1, 6)); - Scatter2DPtr mult = bookScatter2D(ix, 1, 6); + Scatter2DPtr mult; + book(mult, ix, 1, 6); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/GeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } } //@} /// @name Histograms //@{ CounterPtr _nChi1,_nChi2; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2015_I1406939); } diff --git a/analyses/pluginBES/BESIII_2016_I1422780.cc b/analyses/pluginBES/BESIII_2016_I1422780.cc --- a/analyses/pluginBES/BESIII_2016_I1422780.cc +++ b/analyses/pluginBES/BESIII_2016_I1422780.cc @@ -1,185 +1,186 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/Beam.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Jpsi/psi2S baryon decay analysis class BESIII_2016_I1422780 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2016_I1422780); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(Beam(), "Beams"); declare(UnstableParticles(), "UFS"); declare(FinalState(), "FS"); // Book histograms if(fuzzyEquals(sqrtS(),3.1,1e-1)) { - _h_xi = bookHisto1D(2, 1, 1); - _h_sigm = bookHisto1D(2, 1, 2); - _h_sigp = bookHisto1D(2, 1, 3); + book(_h_xi , 2, 1, 1); + book(_h_sigm, 2, 1, 2); + book(_h_sigp, 2, 1, 3); } else if (fuzzyEquals(sqrtS(), 3.686, 1E-1)) { - _h_xi = bookHisto1D(2, 1, 4); - _h_sigm = bookHisto1D(2, 1, 5); - _h_sigp = bookHisto1D(2, 1, 6); + book(_h_xi , 2, 1, 4); + book(_h_sigm, 2, 1, 5); + book(_h_sigp, 2, 1, 6); } } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - nRes[child.pdgId()]-=1; + nRes[child.pid()]-=1; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { // get the axis, direction of incoming electron const ParticlePair& beams = apply(event, "Beams").beams(); Vector3 axis; - if(beams.first.pdgId()>0) + if(beams.first.pid()>0) axis = beams.first .momentum().p3().unit(); else axis = beams.second.momentum().p3().unit(); // types of final state particles const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } const UnstableParticles & ufs = apply(event, "UFS"); - foreach (const Particle& p, ufs.particles(Cuts::abspid==3312 or Cuts::abspid==3224 or Cuts::abspid==3114)) { + for (const Particle& p : ufs.particles(Cuts::abspid==3312 or Cuts::abspid==3224 or Cuts::abspid==3114)) { if(p.children().empty()) continue; map nRes=nCount; int ncount = ntotal; findChildren(p,nRes,ncount); bool matched=false; // check for antiparticle - foreach (const Particle& p2, ufs.particles(Cuts::pid==-p.pdgId())) { + for (const Particle& p2 : ufs.particles(Cuts::pid==-p.pid())) { if(p2.children().empty()) continue; map nRes2=nRes; int ncount2 = ncount; findChildren(p2,nRes2,ncount2); if(ncount2==0) { matched = true; for(auto const & val : nRes2) { if(val.second!=0) { matched = false; break; } } // fond baryon and antibaryon if(matched) { // calc cosine double ctheta; - if(p.pdgId()>0) + if(p.pid()>0) ctheta = p .momentum().p3().unit().dot(axis); else ctheta = p2.momentum().p3().unit().dot(axis); - if(abs(p.pdgId())==3312) - _h_xi ->fill(ctheta,event.weight()); - else if(abs(p.pdgId())==3114) - _h_sigm->fill(ctheta,event.weight()); - else if(abs(p.pdgId())==3224) - _h_sigp->fill(ctheta,event.weight()); + if(abs(p.pid())==3312) + _h_xi ->fill(ctheta); + else if(abs(p.pid())==3114) + _h_sigm->fill(ctheta); + else if(abs(p.pid())==3224) + _h_sigp->fill(ctheta); break; } } } if(matched) break; } } pair > calcAlpha(Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,make_pair(0.,0.)); double d = 3./(pow(hist->xMax(),3)-pow(hist->xMin(),3)); double c = 3.*(hist->xMax()-hist->xMin())/(pow(hist->xMax(),3)-pow(hist->xMin(),3)); double sum1(0.),sum2(0.),sum3(0.),sum4(0.),sum5(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double a = d*(bin.xMax() - bin.xMin()); double b = d/3.*(pow(bin.xMax(),3) - pow(bin.xMin(),3)); double Ei = bin.areaErr(); sum1 += a*Oi/sqr(Ei); sum2 += b*Oi/sqr(Ei); sum3 += sqr(a)/sqr(Ei); sum4 += sqr(b)/sqr(Ei); sum5 += a*b/sqr(Ei); } // calculate alpha double alpha = (-c*sum1 + sqr(c)*sum2 + sum3 - c*sum5)/(sum1 - c*sum2 + c*sum4 - sum5); // and error double cc = -pow((sum3 + sqr(c)*sum4 - 2*c*sum5),3); double bb = -2*sqr(sum3 + sqr(c)*sum4 - 2*c*sum5)*(sum1 - c*sum2 + c*sum4 - sum5); double aa = sqr(sum1 - c*sum2 + c*sum4 - sum5)*(-sum3 - sqr(c)*sum4 + sqr(sum1 - c*sum2 + c*sum4 - sum5) + 2*c*sum5); double dis = sqr(bb)-4.*aa*cc; if(dis>0.) { dis = sqrt(dis); return make_pair(alpha,make_pair(0.5*(-bb+dis)/aa,-0.5*(-bb-dis)/aa)); } else { return make_pair(alpha,make_pair(0.,0.)); } } /// Normalise histograms etc., after the run void finalize() { // find energy int ioff=-1; if(fuzzyEquals(sqrtS(),3.1,1e-1)) { ioff=0; } else if (fuzzyEquals(sqrtS(), 3.686, 1E-1)) { ioff=1; } vector< pair > > alpha; normalize(_h_xi,1.,false); alpha.push_back(calcAlpha(_h_xi)); normalize(_h_sigm,1.,false); alpha.push_back(calcAlpha(_h_sigm)); normalize(_h_sigp,1.,false); alpha.push_back(calcAlpha(_h_sigp)); - Scatter2DPtr _h_alpha = bookScatter2D(1,1,3); + Scatter2DPtr _h_alpha; + book(_h_alpha,1,1,3); for(unsigned int ix=0;ix<3;++ix) _h_alpha->addPoint(1.+ix+3.*ioff, alpha[ix].first, make_pair(0.5,0.5), make_pair(alpha[ix].second.first,alpha[ix].second.second) ); } //@} /// @name Histograms //@{ Histo1DPtr _h_xi,_h_sigm,_h_sigp; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2016_I1422780); } diff --git a/analyses/pluginBES/BESIII_2016_I1495838.cc b/analyses/pluginBES/BESIII_2016_I1495838.cc --- a/analyses/pluginBES/BESIII_2016_I1495838.cc +++ b/analyses/pluginBES/BESIII_2016_I1495838.cc @@ -1,117 +1,117 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2016_I1495838 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2016_I1495838); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); book(_nPsi, "TMP/psi"); } void findChildren(const Particle & p,map & nRes, int &ncount) { for (const Particle &child : p.children()) { if(child.children().empty()) { --nRes[child.pid()]; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); for (const Particle& p: fs.particles()) { nCount[p.pid()] += 1; ++ntotal; } const FinalState& ufs = apply(event, "UFS"); // find the psis - foreach (const Particle& p, ufs.particles(Cuts::pid==443)) { + for (const Particle& p : ufs.particles(Cuts::pid==443)) { if(p.children().empty()) continue; map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); // psi pi+pi- if(ncount!=2) continue; bool matched = true; for(auto const & val : nRes) { if(abs(val.first)==211) { if(val.second !=1) { matched = false; break; } } else if(val.second!=0) { matched = false; break; } } if(matched) { _nPsi->fill(); break; } } } /// Normalise histograms etc., after the run void finalize() { double fact = crossSection()/ sumOfWeights() /picobarn; double sigma = _nPsi->val()*fact; double error = _nPsi->err()*fact; for(unsigned int ix=1;ix<3;++ix) { Scatter2D temphisto(refData(ix, 1, 1)); Scatter2DPtr mult; book(mult, ix, 1, 1); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/GeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } } //@} /// @name Histograms //@{ CounterPtr _nPsi; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2016_I1495838); } diff --git a/analyses/pluginBES/BESIII_2017_I1506414.cc b/analyses/pluginBES/BESIII_2017_I1506414.cc --- a/analyses/pluginBES/BESIII_2017_I1506414.cc +++ b/analyses/pluginBES/BESIII_2017_I1506414.cc @@ -1,181 +1,183 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/Beam.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2017_I1506414 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2017_I1506414); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(Beam(), "Beams"); declare(UnstableParticles(), "UFS"); declare(FinalState(), "FS"); // Book histograms if(fuzzyEquals(sqrtS(),3.1,1e-1)) { - _h_xi = bookHisto1D(1, 1, 2); - _h_sig = bookHisto1D(1, 1, 1); + book(_h_xi , 1, 1, 2); + book(_h_sig, 1, 1, 1); } else if (fuzzyEquals(sqrtS(), 3.686, 1E-1)) { - _h_xi = bookHisto1D(1, 1, 4); - _h_sig = bookHisto1D(1, 1, 3); + book(_h_xi , 1, 1, 4); + book(_h_sig, 1, 1, 3); } } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - nRes[child.pdgId()]-=1; + nRes[child.pid()]-=1; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { // get the axis, direction of incoming electron const ParticlePair& beams = apply(event, "Beams").beams(); Vector3 axis; - if(beams.first.pdgId()>0) + if(beams.first.pid()>0) axis = beams.first .momentum().p3().unit(); else axis = beams.second.momentum().p3().unit(); // types of final state particles const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } const UnstableParticles & ufs = apply(event, "UFS"); - foreach (const Particle& p, ufs.particles(Cuts::abspid==3322 or Cuts::abspid==3214)) { + for (const Particle& p : ufs.particles(Cuts::abspid==3322 or Cuts::abspid==3214)) { if(p.children().empty()) continue; map nRes=nCount; int ncount = ntotal; findChildren(p,nRes,ncount); bool matched=false; // check for antiparticle - foreach (const Particle& p2, ufs.particles(Cuts::pid==-p.pdgId())) { + for (const Particle& p2 : ufs.particles(Cuts::pid==-p.pid())) { if(p2.children().empty()) continue; map nRes2=nRes; int ncount2 = ncount; findChildren(p2,nRes2,ncount2); if(ncount2==0) { matched = true; for(auto const & val : nRes2) { if(val.second!=0) { matched = false; break; } } // fond baryon and antibaryon if(matched) { // calc cosine double ctheta; - if(p.pdgId()>0) + if(p.pid()>0) ctheta = p .momentum().p3().unit().dot(axis); else ctheta = p2.momentum().p3().unit().dot(axis); - if(abs(p.pdgId())==3322) - _h_xi ->fill(ctheta,event.weight()); - else if(abs(p.pdgId())==3214) - _h_sig->fill(ctheta,event.weight()); + if(abs(p.pid())==3322) + _h_xi ->fill(ctheta); + else if(abs(p.pid())==3214) + _h_sig->fill(ctheta); break; } } } if(matched) break; } } pair > calcAlpha(Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,make_pair(0.,0.)); double d = 3./(pow(hist->xMax(),3)-pow(hist->xMin(),3)); double c = 3.*(hist->xMax()-hist->xMin())/(pow(hist->xMax(),3)-pow(hist->xMin(),3)); double sum1(0.),sum2(0.),sum3(0.),sum4(0.),sum5(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double a = d*(bin.xMax() - bin.xMin()); double b = d/3.*(pow(bin.xMax(),3) - pow(bin.xMin(),3)); double Ei = bin.areaErr(); sum1 += a*Oi/sqr(Ei); sum2 += b*Oi/sqr(Ei); sum3 += sqr(a)/sqr(Ei); sum4 += sqr(b)/sqr(Ei); sum5 += a*b/sqr(Ei); } // calculate alpha double alpha = (-c*sum1 + sqr(c)*sum2 + sum3 - c*sum5)/(sum1 - c*sum2 + c*sum4 - sum5); // and error double cc = -pow((sum3 + sqr(c)*sum4 - 2*c*sum5),3); double bb = -2*sqr(sum3 + sqr(c)*sum4 - 2*c*sum5)*(sum1 - c*sum2 + c*sum4 - sum5); double aa = sqr(sum1 - c*sum2 + c*sum4 - sum5)*(-sum3 - sqr(c)*sum4 + sqr(sum1 - c*sum2 + c*sum4 - sum5) + 2*c*sum5); double dis = sqr(bb)-4.*aa*cc; if(dis>0.) { dis = sqrt(dis); return make_pair(alpha,make_pair(0.5*(-bb+dis)/aa,-0.5*(-bb-dis)/aa)); } else { return make_pair(alpha,make_pair(0.,0.)); } } /// Normalise histograms etc., after the run void finalize() { // find energy int ioff=-1; if(fuzzyEquals(sqrtS(),3.1,1e-1)) { ioff=0; } else if (fuzzyEquals(sqrtS(), 3.686, 1E-1)) { ioff=1; } normalize(_h_xi,1.,false); - Scatter2DPtr _h_alpha_xi = bookScatter2D(2,2*ioff+2,1); + Scatter2DPtr _h_alpha_xi; + book(_h_alpha_xi,2,2*ioff+2,1); pair > alpha = calcAlpha(_h_xi); _h_alpha_xi->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second.first,alpha.second.second) ); normalize(_h_sig,1.,false); - Scatter2DPtr _h_alpha_sig = bookScatter2D(2,2*ioff+1,1); + Scatter2DPtr _h_alpha_sig; + book(_h_alpha_sig,2,2*ioff+1,1); alpha = calcAlpha(_h_sig); _h_alpha_sig->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second.first,alpha.second.second) ); } //@} /// @name Histograms //@{ Histo1DPtr _h_xi,_h_sig; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2017_I1506414); } diff --git a/analyses/pluginBES/BESIII_2017_I1510563.cc b/analyses/pluginBES/BESIII_2017_I1510563.cc --- a/analyses/pluginBES/BESIII_2017_I1510563.cc +++ b/analyses/pluginBES/BESIII_2017_I1510563.cc @@ -1,176 +1,180 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/Beam.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Jpsi/psi2S baryon decay analysis class BESIII_2017_I1510563 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2017_I1510563); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(Beam(), "Beams"); declare(UnstableParticles(), "UFS"); declare(FinalState(), "FS"); // Book histograms if(fuzzyEquals(sqrtS(),3.1,1e-1)) { - _h_lam = bookHisto1D(1, 1, 1); - _h_sig = bookHisto1D(1, 1, 3); + book(_h_lam, 1, 1, 1); + book(_h_sig, 1, 1, 3); } else if (fuzzyEquals(sqrtS(), 3.686, 1E-1)) { - _h_lam = bookHisto1D(1, 1, 2); - _h_sig = bookHisto1D(1, 1, 4); + book(_h_lam, 1, 1, 2); + book(_h_sig, 1, 1, 4); } } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - nRes[child.pdgId()]-=1; + nRes[child.pid()]-=1; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { // get the axis, direction of incoming electron const ParticlePair& beams = apply(event, "Beams").beams(); Vector3 axis; - if(beams.first.pdgId()>0) + if(beams.first.pid()>0) axis = beams.first .momentum().p3().unit(); else axis = beams.second.momentum().p3().unit(); // types of final state particles const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } // loop over lambda0 and sigma0 baryons const UnstableParticles & ufs = apply(event, "UFS"); - foreach (const Particle& p, ufs.particles(Cuts::abspid==3122 or Cuts::abspid==3212)) { + for (const Particle& p : ufs.particles(Cuts::abspid==3122 or Cuts::abspid==3212)) { if(p.children().empty()) continue; map nRes=nCount; int ncount = ntotal; findChildren(p,nRes,ncount); bool matched=false; // check for antiparticle - foreach (const Particle& p2, ufs.particles(Cuts::pid==-p.pdgId())) { + for (const Particle& p2 : ufs.particles(Cuts::pid==-p.pid())) { if(p2.children().empty()) continue; map nRes2=nRes; int ncount2 = ncount; findChildren(p2,nRes2,ncount2); if(ncount2==0) { matched = true; for(auto const & val : nRes2) { if(val.second!=0) { matched = false; break; } } // fond baryon and antibaryon if(matched) { // calc cosine double ctheta; - if(p.pdgId()>0) + if(p.pid()>0) ctheta = p .momentum().p3().unit().dot(axis); else ctheta = p2.momentum().p3().unit().dot(axis); - if(abs(p.pdgId())==3122) - _h_lam->fill(ctheta,event.weight()); + if(abs(p.pid())==3122) + _h_lam->fill(ctheta); else - _h_sig->fill(ctheta,event.weight()); + _h_sig->fill(ctheta); break; } } } if(matched) break; } } pair > calcAlpha(Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,make_pair(0.,0.)); double sum1(0.),sum2(0.),sum3(0.),sum4(0.),sum5(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double a = 1.5*(bin.xMax() - bin.xMin()); double b = 0.5*(pow(bin.xMax(),3) - pow(bin.xMin(),3)); double Ei = bin.areaErr(); sum1 += a*Oi/sqr(Ei); sum2 += b*Oi/sqr(Ei); sum3 += sqr(a)/sqr(Ei); sum4 += sqr(b)/sqr(Ei); sum5 += a*b/sqr(Ei); } // calculate alpha double alpha = (-3*sum1 + 9*sum2 + sum3 - 3*sum5)/(sum1 - 3*sum2 + 3*sum4 - sum5); // and error double cc = -pow((sum3 + 9*sum4 - 6*sum5),3); double bb = -2*sqr(sum3 + 9*sum4 - 6*sum5)*(sum1 - 3*sum2 + 3*sum4 - sum5); double aa = sqr(sum1 - 3*sum2 + 3*sum4 - sum5)*(-sum3 - 9*sum4 + sqr(sum1 - 3*sum2 + 3*sum4 - sum5) + 6*sum5); double dis = sqr(bb)-4.*aa*cc; if(dis>0.) { dis = sqrt(dis); return make_pair(alpha,make_pair(0.5*(-bb+dis)/aa,-0.5*(-bb-dis)/aa)); } else { return make_pair(alpha,make_pair(0.,0.)); } } /// Normalise histograms etc., after the run void finalize() { // find energy int ioff=-1; if(fuzzyEquals(sqrtS(),3.1,1e-1)) { ioff=0; } else if (fuzzyEquals(sqrtS(), 3.686, 1E-1)) { ioff=1; } // normalize normalize(_h_lam); pair > alpha = calcAlpha(_h_lam); - Scatter2DPtr _h_alpha_lam = bookScatter2D(2,2*ioff+1,1); - _h_alpha_lam->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second.first,alpha.second.second) ); + Scatter2DPtr _h_alpha_lam; + book(_h_alpha_lam, 2,2*ioff+1,1); + _h_alpha_lam->addPoint(0.5, alpha.first, make_pair(0.5,0.5), + make_pair(alpha.second.first,alpha.second.second) ); normalize(_h_sig); alpha = calcAlpha(_h_sig); - Scatter2DPtr _h_alpha_sig = bookScatter2D(2,2*ioff+2,1); - _h_alpha_sig->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second.first,alpha.second.second) ); + Scatter2DPtr _h_alpha_sig; + book(_h_alpha_sig, 2,2*ioff+2,1); + _h_alpha_sig->addPoint(0.5, alpha.first, make_pair(0.5,0.5), + make_pair(alpha.second.first,alpha.second.second) ); } //@} /// @name Histograms //@{ Histo1DPtr _h_lam,_h_sig; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2017_I1510563); } diff --git a/analyses/pluginBES/BESIII_2017_I1519425.cc b/analyses/pluginBES/BESIII_2017_I1519425.cc --- a/analyses/pluginBES/BESIII_2017_I1519425.cc +++ b/analyses/pluginBES/BESIII_2017_I1519425.cc @@ -1,98 +1,98 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2017_I1519425 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2017_I1519425); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS"); // Book histograms - _h_q2_K = bookHisto1D(1, 1, 1); - _h_q2_pi = bookHisto1D(2, 1, 1); + book(_h_q2_K , 1, 1, 1); + book(_h_q2_pi, 2, 1, 1); } // Calculate the Q2 using mother and daugher meson double q2(const Particle& B, int mesonID) { FourMomentum q = B.mom() - filter_select(B.children(), Cuts::pid==mesonID)[0]; return q*q; } // Check for explicit decay into pdgids bool isSemileptonicDecay(const Particle& mother, vector ids) { // Trivial check to ignore any other decays but the one in question modulo photons const Particles children = mother.children(Cuts::pid!=PID::PHOTON); if (children.size()!=ids.size()) return false; // Check for the explicit decay return all(ids, [&](int i){return count(children, hasPID(i))==1;}); } /// Perform the per-event analysis void analyze(const Event& event) { // Loop over D+/- mesons - foreach(const Particle& p, apply(event, "UFS").particles(Cuts::abspid==PID::DPLUS)) { + for (const Particle& p : apply(event, "UFS").particles(Cuts::abspid==PID::DPLUS)) { if (isSemileptonicDecay(p, {PID::PI0, PID::POSITRON, PID::NU_E}) || isSemileptonicDecay(p, {PID::PI0, PID::ELECTRON, PID::NU_EBAR})) { - _h_q2_pi->fill(q2(p, PID::PI0), event.weight()); + _h_q2_pi->fill(q2(p, PID::PI0)); } else if(isSemileptonicDecay(p, {-311, PID::POSITRON, PID::NU_E})) { - _h_q2_K ->fill(q2(p, -311), event.weight()); + _h_q2_K ->fill(q2(p, -311)); } else if(isSemileptonicDecay(p, { 311, PID::ELECTRON, PID::NU_EBAR})) { - _h_q2_K ->fill(q2(p, 311), event.weight()); + _h_q2_K ->fill(q2(p, 311)); } else if(isSemileptonicDecay(p, {PID::K0S, PID::POSITRON, PID::NU_E}) || isSemileptonicDecay(p, {PID::K0S, PID::ELECTRON, PID::NU_EBAR})) { - _h_q2_K ->fill(q2(p, PID::K0S), event.weight()); + _h_q2_K ->fill(q2(p, PID::K0S)); } else if(isSemileptonicDecay(p, {PID::K0L, PID::POSITRON, PID::NU_E}) || isSemileptonicDecay(p, {PID::K0L, PID::ELECTRON, PID::NU_EBAR})) { - _h_q2_K ->fill(q2(p, PID::K0L), event.weight()); + _h_q2_K ->fill(q2(p, PID::K0L)); } } } /// Normalise histograms etc., after the run void finalize() { // normalize to unity normalize(_h_q2_K ); normalize(_h_q2_pi); } //@} /// @name Histograms //@{ Histo1DPtr _h_q2_K,_h_q2_pi; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2017_I1519425); } diff --git a/analyses/pluginBES/BESIII_2017_I1628093.cc b/analyses/pluginBES/BESIII_2017_I1628093.cc --- a/analyses/pluginBES/BESIII_2017_I1628093.cc +++ b/analyses/pluginBES/BESIII_2017_I1628093.cc @@ -1,126 +1,127 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2017_I1628093 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2017_I1628093); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); - _nLambda = bookCounter("/TMP/nLambda"); + book(_nLambda,"/TMP/nLambda"); } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - nRes[child.pdgId()]-=1; + nRes[child.pid()]-=1; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } const FinalState& ufs = apply(event, "UFS"); for(unsigned int ix=0;ix nRes = nCount; int ncount = ntotal; findChildren(p1,nRes,ncount); bool matched=false; for(unsigned int iy=0;iy nRes2 = nRes; int ncount2 = ncount; findChildren(p2,nRes2,ncount2); if(ncount2!=0) continue; matched=true; for(auto const & val : nRes2) { if(val.second!=0) { matched = false; break; } } if(matched) { break; } } if(matched) { - _nLambda->fill(event.weight()); + _nLambda->fill(); break; } } } /// Normalise histograms etc., after the run void finalize() { double sigma = _nLambda->val(); double error = _nLambda->err(); sigma *= crossSection()/ sumOfWeights() /picobarn; error *= crossSection()/ sumOfWeights() /picobarn; Scatter2D temphisto(refData(1, 1, 1)); - Scatter2DPtr mult = bookScatter2D(1, 1, 1); + Scatter2DPtr mult; + book(mult,1, 1, 1); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/MeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } //@} /// @name Histograms //@{ CounterPtr _nLambda; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2017_I1628093); } diff --git a/analyses/pluginBES/BESIII_2018_I1641075.cc b/analyses/pluginBES/BESIII_2018_I1641075.cc --- a/analyses/pluginBES/BESIII_2018_I1641075.cc +++ b/analyses/pluginBES/BESIII_2018_I1641075.cc @@ -1,94 +1,94 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2018_I1641075 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2018_I1641075); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS"); // Book histograms - _h_m = bookHisto1D(1, 1, 5); + book(_h_m, 1, 1, 5); } void findDecayProducts(const Particle & mother, unsigned int & nstable, unsigned int & ngamma, unsigned int & npip, unsigned int & npim, FourMomentum & ptot) { for(const Particle & p : mother.children()) { - int id = p.pdgId(); + int id = p.pid(); if (id == PID::PIMINUS ) { ++npim; ++nstable; ptot += p.momentum(); } else if (id == PID::PIPLUS) { ++npip; ++nstable; ptot += p.momentum(); } else if ( !p.children().empty() ) { findDecayProducts(p, nstable, ngamma,npip,npim,ptot); } else if (id == PID::GAMMA) { ++ngamma; ++nstable; } else ++nstable; } } /// Perform the per-event analysis void analyze(const Event& event) { // Loop over eta' mesons - foreach(const Particle& p, apply(event, "UFS").particles(Cuts::pid==331)) { + for (const Particle& p : apply(event, "UFS").particles(Cuts::pid==331)) { unsigned nstable(0),ngamma(0),npip(0),npim(0); FourMomentum ptot; findDecayProducts(p,nstable,ngamma,npip,npim,ptot); if(nstable==3 && npim==1 && npip==1 && ngamma==1) - _h_m->fill(ptot.mass(),event.weight()); + _h_m->fill(ptot.mass()); } } /// Normalise histograms etc., after the run void finalize() { normalize(_h_m); } //@} /// @name Histograms //@{ Histo1DPtr _h_m; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2018_I1641075); } diff --git a/analyses/pluginBES/BESIII_2018_I1658762.cc b/analyses/pluginBES/BESIII_2018_I1658762.cc --- a/analyses/pluginBES/BESIII_2018_I1658762.cc +++ b/analyses/pluginBES/BESIII_2018_I1658762.cc @@ -1,127 +1,131 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/Beam.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief psi2S baryon decay analysis class BESIII_2018_I1658762 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2018_I1658762); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(Beam(), "Beams"); declare(UnstableParticles(), "UFS"); declare(FinalState(), "FS"); // Book histograms - _h_proton = bookHisto1D("ctheta_p",20,-1.,1.); - _h_neutron = bookHisto1D("ctheta_n",20,-1.,1.); + book(_h_proton ,"ctheta_p",20,-1.,1.); + book(_h_neutron,"ctheta_n",20,-1.,1.); } /// Perform the per-event analysis void analyze(const Event& event) { // get the axis, direction of incoming electron const ParticlePair& beams = apply(event, "Beams").beams(); Vector3 axis; - if(beams.first.pdgId()>0) + if(beams.first.pid()>0) axis = beams.first .momentum().p3().unit(); else axis = beams.second.momentum().p3().unit(); // types of final state particles const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); Particle outgoing; - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; - if(p.pdgId()==2212 || p.pdgId()==2112) + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; + if(p.pid()==2212 || p.pid()==2112) outgoing = p; ++ntotal; } if(ntotal==2) { if(nCount[2212]==1 && nCount[-2212]==1) { - _h_proton->fill(outgoing.momentum().p3().unit().dot(axis),event.weight()); + _h_proton->fill(outgoing.momentum().p3().unit().dot(axis)); } else if(nCount[2112]==1 && nCount[-2112]==1) { - _h_neutron->fill(outgoing.momentum().p3().unit().dot(axis),event.weight()); + _h_neutron->fill(outgoing.momentum().p3().unit().dot(axis)); } } } pair > calcAlpha(Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,make_pair(0.,0.)); double sum1(0.),sum2(0.),sum3(0.),sum4(0.),sum5(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double a = 1.5*(bin.xMax() - bin.xMin()); double b = 0.5*(pow(bin.xMax(),3) - pow(bin.xMin(),3)); double Ei = bin.areaErr(); sum1 += a*Oi/sqr(Ei); sum2 += b*Oi/sqr(Ei); sum3 += sqr(a)/sqr(Ei); sum4 += sqr(b)/sqr(Ei); sum5 += a*b/sqr(Ei); } // calculate alpha double alpha = (-3*sum1 + 9*sum2 + sum3 - 3*sum5)/(sum1 - 3*sum2 + 3*sum4 - sum5); // and error double cc = -pow((sum3 + 9*sum4 - 6*sum5),3); double bb = -2*sqr(sum3 + 9*sum4 - 6*sum5)*(sum1 - 3*sum2 + 3*sum4 - sum5); double aa = sqr(sum1 - 3*sum2 + 3*sum4 - sum5)*(-sum3 - 9*sum4 + sqr(sum1 - 3*sum2 + 3*sum4 - sum5) + 6*sum5); double dis = sqr(bb)-4.*aa*cc; if(dis>0.) { dis = sqrt(dis); return make_pair(alpha,make_pair(0.5*(-bb+dis)/aa,-0.5*(-bb-dis)/aa)); } else { return make_pair(alpha,make_pair(0.,0.)); } } /// Normalise histograms etc., after the run void finalize() { // proton normalize(_h_proton ); pair > alpha = calcAlpha(_h_proton); - Scatter2DPtr _h_alpha_proton = bookScatter2D(1,1,1); - _h_alpha_proton->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second.first,alpha.second.second) ); + Scatter2DPtr _h_alpha_proton; + book(_h_alpha_proton,1,1,1); + _h_alpha_proton->addPoint(0.5, alpha.first, make_pair(0.5,0.5), + make_pair(alpha.second.first,alpha.second.second) ); // neutron normalize(_h_neutron); alpha = calcAlpha(_h_neutron); - Scatter2DPtr _h_alpha_neutron = bookScatter2D(1,1,2); - _h_alpha_neutron->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second.first,alpha.second.second) ); + Scatter2DPtr _h_alpha_neutron; + book(_h_alpha_neutron,1,1,2); + _h_alpha_neutron->addPoint(0.5, alpha.first, make_pair(0.5,0.5), + make_pair(alpha.second.first,alpha.second.second) ); } //@} /// @name Histograms //@{ Histo1DPtr _h_proton,_h_neutron; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2018_I1658762); } diff --git a/analyses/pluginBES/BESIII_2018_I1697377.cc b/analyses/pluginBES/BESIII_2018_I1697377.cc --- a/analyses/pluginBES/BESIII_2018_I1697377.cc +++ b/analyses/pluginBES/BESIII_2018_I1697377.cc @@ -1,92 +1,92 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2018_I1697377 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2018_I1697377); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS"); // Book histograms - _h_m = bookHisto1D(1, 1, 1); + book(_h_m, 1, 1, 1); } void findDecayProducts(const Particle & mother, unsigned int & nstable, unsigned int & neta, unsigned int & nep, unsigned int & nem, FourMomentum & ptot) { for(const Particle & p : mother.children()) { - int id = p.pdgId(); + int id = p.pid(); if (id == PID::EMINUS ) { ++nem; ++nstable; ptot += p.momentum(); } else if (id == PID::EPLUS) { ++nep; ++nstable; ptot += p.momentum(); } else if (id == PID::ETA) { ++neta; ++nstable; } else if ( !p.children().empty() ) { findDecayProducts(p, nstable, neta,nep,nem,ptot); } else ++nstable; } } /// Perform the per-event analysis void analyze(const Event& event) { // Loop over J/psi mesons - foreach(const Particle& p, apply(event, "UFS").particles(Cuts::pid==443)) { + for (const Particle& p : apply(event, "UFS").particles(Cuts::pid==443)) { unsigned nstable(0),neta(0),nep(0),nem(0); FourMomentum ptot; findDecayProducts(p,nstable,neta,nep,nem,ptot); if(nstable==3 && nem==1 && nem==1 && neta==1) { - _h_m->fill(ptot.mass(),event.weight()); + _h_m->fill(ptot.mass()); } } } /// Normalise histograms etc., after the run void finalize() { normalize(_h_m); // normalize to unity } //@} /// @name Histograms //@{ Histo1DPtr _h_m; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2018_I1697377); } diff --git a/analyses/pluginBES/BESIII_2019_I1691850.cc b/analyses/pluginBES/BESIII_2019_I1691850.cc --- a/analyses/pluginBES/BESIII_2019_I1691850.cc +++ b/analyses/pluginBES/BESIII_2019_I1691850.cc @@ -1,381 +1,386 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/Beam.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief JPsi > Lambda, Lambdabar class BESIII_2019_I1691850 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2019_I1691850); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(Beam(), "Beams"); declare(UnstableParticles(), "UFS"); declare(FinalState(), "FS"); // Book histograms - _h_T1_p = bookHisto1D("T1_p",20,-1.,1.); - _h_T2_p = bookHisto1D("T2_p",20,-1.,1.); - _h_T3_p = bookHisto1D("T3_p",20,-1.,1.); - _h_T4_p = bookHisto1D("T4_p",20,-1.,1.); - _h_T5_p = bookHisto1D("T5_p",20,-1.,1.); - - _h_T1_n = bookHisto1D("T1_n",20,-1.,1.); - _h_T2_n = bookHisto1D("T2_n",20,-1.,1.); - _h_T3_n = bookHisto1D("T3_n",20,-1.,1.); - _h_T4_n = bookHisto1D("T4_n",20,-1.,1.); - _h_T5_n = bookHisto1D("T5_n",20,-1.,1.); + book(_h_T1_p, "T1_p",20,-1.,1.); + book(_h_T2_p, "T2_p",20,-1.,1.); + book(_h_T3_p, "T3_p",20,-1.,1.); + book(_h_T4_p, "T4_p",20,-1.,1.); + book(_h_T5_p, "T5_p",20,-1.,1.); + + book(_h_T1_n, "T1_n",20,-1.,1.); + book(_h_T2_n, "T2_n",20,-1.,1.); + book(_h_T3_n, "T3_n",20,-1.,1.); + book(_h_T4_n, "T4_n",20,-1.,1.); + book(_h_T5_n, "T5_n",20,-1.,1.); - _h_cThetaL = bookHisto1D("cThetaL",20,-1.,1.); + book(_h_cThetaL,"cThetaL",20,-1.,1.); - _h_mu_p = bookHisto1D("mu_p",20,-1.,1.); - _h_mu_n = bookHisto1D("mu_n",20,-1.,1.); - _wsum_p = 0.; - _wsum_n = 0.; + book(_h_mu_p, "mu_p",20,-1.,1.); + book(_h_mu_n, "mu_n",20,-1.,1.); + book(_wsum_p,"TMP/wsum_p"); + book(_wsum_n,"TMP/wsum_n"); } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - nRes[child.pdgId()]-=1; + nRes[child.pid()]-=1; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { // get the axis, direction of incoming electron const ParticlePair& beams = apply(event, "Beams").beams(); Vector3 axis; - if(beams.first.pdgId()>0) + if(beams.first.pid()>0) axis = beams.first .momentum().p3().unit(); else axis = beams.second.momentum().p3().unit(); // types of final state particles const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } // loop over lambda0 baryons const UnstableParticles & ufs = apply(event, "UFS"); Particle Lambda,LamBar; bool matched(false); - foreach (const Particle& p, ufs.particles(Cuts::abspid==3122)) { + for (const Particle& p : ufs.particles(Cuts::abspid==3122)) { if(p.children().empty()) continue; map nRes=nCount; int ncount = ntotal; findChildren(p,nRes,ncount); matched=false; // check for antiparticle - foreach (const Particle& p2, ufs.particles(Cuts::pid==-p.pdgId())) { + for (const Particle& p2 : ufs.particles(Cuts::pid==-p.pid())) { if(p2.children().empty()) continue; map nRes2=nRes; int ncount2 = ncount; findChildren(p2,nRes2,ncount2); if(ncount2==0) { matched = true; for(auto const & val : nRes2) { if(val.second!=0) { matched = false; break; } } // fond baryon and antibaryon if(matched) { - if(p.pdgId()>0) { + if(p.pid()>0) { Lambda = p; LamBar = p2; } else { Lambda = p2; LamBar = p; } break; } } } if(matched) break; } if(!matched) vetoEvent; Particle proton; matched = false; for (const Particle & p : Lambda.children()) { - if(p.pdgId()==2212) { + if(p.pid()==2212) { matched=true; proton=p; } } if(!matched) vetoEvent; Particle baryon; int mode(-1); for (const Particle & p : LamBar.children()) { - if(p.pdgId()==-2212) { + if(p.pid()==-2212) { baryon=p; mode=0; } - else if(p.pdgId()==-2112) { + else if(p.pid()==-2112) { baryon=p; mode=1; } } // boost to the Lambda rest frame LorentzTransform boost1 = LorentzTransform::mkFrameTransformFromBeta(Lambda.momentum().betaVec()); Vector3 e1z = Lambda.momentum().p3().unit(); Vector3 e1y = e1z.cross(axis).unit(); Vector3 e1x = e1y.cross(e1z).unit(); Vector3 axis1 = boost1.transform(proton.momentum()).p3().unit(); double n1x(e1x.dot(axis1)),n1y(e1y.dot(axis1)),n1z(e1z.dot(axis1)); // boost to the Lambda bar LorentzTransform boost2 = LorentzTransform::mkFrameTransformFromBeta(LamBar.momentum().betaVec()); Vector3 axis2 = boost2.transform(baryon.momentum()).p3().unit(); double n2x(e1x.dot(axis2)),n2y(e1y.dot(axis2)),n2z(e1z.dot(axis2)); double cosL = axis.dot(Lambda.momentum().p3().unit()); double sinL = sqrt(1.-sqr(cosL)); double T1 = sqr(sinL)*n1x*n2x+sqr(cosL)*n1z*n2z; double T2 = -sinL*cosL*(n1x*n2z+n1z*n2x); double T3 = -sinL*cosL*n1y; double T4 = -sinL*cosL*n2y; double T5 = n1z*n2z-sqr(sinL)*n1y*n2y; double mu = n1y-n2y; if(mode==0) { - _h_T1_p->fill(cosL,T1*event.weight()); - _h_T2_p->fill(cosL,T2*event.weight()); - _h_T3_p->fill(cosL,T3*event.weight()); - _h_T4_p->fill(cosL,T4*event.weight()); - _h_T5_p->fill(cosL,T5*event.weight()); - _h_mu_p->fill(cosL,mu*event.weight()); - _wsum_p+=event.weight(); + _h_T1_p->fill(cosL,T1); + _h_T2_p->fill(cosL,T2); + _h_T3_p->fill(cosL,T3); + _h_T4_p->fill(cosL,T4); + _h_T5_p->fill(cosL,T5); + _h_mu_p->fill(cosL,mu); + _wsum_p->fill(); } else { - _h_T1_n->fill(cosL,T1*event.weight()); - _h_T2_n->fill(cosL,T2*event.weight()); - _h_T3_n->fill(cosL,T3*event.weight()); - _h_T4_n->fill(cosL,T4*event.weight()); - _h_T5_n->fill(cosL,T5*event.weight()); - _h_mu_n->fill(cosL,mu*event.weight()); - _wsum_n+=event.weight(); + _h_T1_n->fill(cosL,T1); + _h_T2_n->fill(cosL,T2); + _h_T3_n->fill(cosL,T3); + _h_T4_n->fill(cosL,T4); + _h_T5_n->fill(cosL,T5); + _h_mu_n->fill(cosL,mu); + _wsum_n->fill(); } - _h_cThetaL->fill(cosL,event.weight()); + _h_cThetaL->fill(cosL); } pair > calcAlpha0(Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,make_pair(0.,0.)); double d = 3./(pow(hist->xMax(),3)-pow(hist->xMin(),3)); double c = 3.*(hist->xMax()-hist->xMin())/(pow(hist->xMax(),3)-pow(hist->xMin(),3)); double sum1(0.),sum2(0.),sum3(0.),sum4(0.),sum5(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double a = d*(bin.xMax() - bin.xMin()); double b = d/3.*(pow(bin.xMax(),3) - pow(bin.xMin(),3)); double Ei = bin.areaErr(); sum1 += a*Oi/sqr(Ei); sum2 += b*Oi/sqr(Ei); sum3 += sqr(a)/sqr(Ei); sum4 += sqr(b)/sqr(Ei); sum5 += a*b/sqr(Ei); } // calculate alpha double alpha = (-c*sum1 + sqr(c)*sum2 + sum3 - c*sum5)/(sum1 - c*sum2 + c*sum4 - sum5); // and error double cc = -pow((sum3 + sqr(c)*sum4 - 2*c*sum5),3); double bb = -2*sqr(sum3 + sqr(c)*sum4 - 2*c*sum5)*(sum1 - c*sum2 + c*sum4 - sum5); double aa = sqr(sum1 - c*sum2 + c*sum4 - sum5)*(-sum3 - sqr(c)*sum4 + sqr(sum1 - c*sum2 + c*sum4 - sum5) + 2*c*sum5); double dis = sqr(bb)-4.*aa*cc; if(dis>0.) { dis = sqrt(dis); return make_pair(alpha,make_pair(0.5*(-bb+dis)/aa,-0.5*(-bb-dis)/aa)); } else { return make_pair(alpha,make_pair(0.,0.)); } } pair calcCoeff(unsigned int imode,Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,0.); double sum1(0.),sum2(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double ai(0.),bi(0.); if(imode==0) { bi = (pow(1.-sqr(bin.xMin()),1.5) - pow(1.-sqr(bin.xMax()),1.5))/3.; } else if(imode>=2 && imode<=4) { bi = ( pow(bin.xMin(),3)*( -5. + 3.*sqr(bin.xMin())) + pow(bin.xMax(),3)*( 5. - 3.*sqr(bin.xMax())))/15.; } else assert(false); double Ei = bin.areaErr(); sum1 += sqr(bi/Ei); sum2 += bi/sqr(Ei)*(Oi-ai); } return make_pair(sum2/sum1,sqrt(1./sum1)); } /// Normalise histograms etc., after the run void finalize() { normalize(_h_cThetaL); - scale(_h_T1_p,1./_wsum_p); - scale(_h_T2_p,1./_wsum_p); - scale(_h_T3_p,1./_wsum_p); - scale(_h_T4_p,1./_wsum_p); - scale(_h_T5_p,1./_wsum_p); + scale(_h_T1_p,1./ *_wsum_p); + scale(_h_T2_p,1./ *_wsum_p); + scale(_h_T3_p,1./ *_wsum_p); + scale(_h_T4_p,1./ *_wsum_p); + scale(_h_T5_p,1./ *_wsum_p); - scale(_h_T1_n,1./_wsum_n); - scale(_h_T2_n,1./_wsum_n); - scale(_h_T3_n,1./_wsum_n); - scale(_h_T4_n,1./_wsum_n); - scale(_h_T5_n,1./_wsum_n); + scale(_h_T1_n,1./ *_wsum_n); + scale(_h_T2_n,1./ *_wsum_n); + scale(_h_T3_n,1./ *_wsum_n); + scale(_h_T4_n,1./ *_wsum_n); + scale(_h_T5_n,1./ *_wsum_n); - scale(_h_mu_p,1./_wsum_p); - scale(_h_mu_n,1./_wsum_n); + scale(_h_mu_p,1./ *_wsum_p); + scale(_h_mu_n,1./ *_wsum_n); // calculate alpha0 pair > alpha0 = calcAlpha0(_h_cThetaL); - Scatter2DPtr _h_alpha0 = bookScatter2D(1,1,1); + Scatter2DPtr _h_alpha0; + book(_h_alpha0,1,1,1); _h_alpha0->addPoint(0.5, alpha0.first, make_pair(0.5,0.5), make_pair(alpha0.second.first,alpha0.second.second) ); double s2 = -1. + sqr(alpha0.first); double s3 = 3 + alpha0.first; double s1 = sqr(s3); // alpha- and alpha+ from proton data pair c_T2_p = calcCoeff(2,_h_T2_p); pair c_T3_p = calcCoeff(3,_h_T3_p); pair c_T4_p = calcCoeff(4,_h_T4_p); double s4 = sqr(c_T2_p.first); double s5 = sqr(c_T3_p.first); double s6 = sqr(c_T4_p.first); double disc = s1*s5*s6*(-9.*s2*s4 + 4.*s1*s5*s6); if(disc>=0.) { disc = sqrt(disc); double aM = sqrt(-1./s2/s6*(2.*s1*s5*s6+disc)); double aP = c_T4_p.first/c_T3_p.first*aM; double aM_M = (2*(alpha0.first*c_T4_p.first*alpha0.second.first + c_T4_p.second*s2)*(disc + 2*s1*s5*s6) - c_T4_p.first*s2*(4*s3*c_T3_p.first*c_T4_p.first*(c_T3_p.first*c_T4_p.first*alpha0.second.first +s3*c_T4_p.first*c_T3_p.second +s3*c_T3_p.first*c_T4_p.second) + (disc*(- 9*s2*s3*c_T2_p.first*c_T3_p.first*c_T4_p.first* c_T2_p.second + 9*((1 - alpha0.first*(3 + 2*alpha0.first))* c_T3_p.first*c_T4_p.first*alpha0.second.first - s2*s3*c_T4_p.first*c_T3_p.second - s2*s3*c_T3_p.first*c_T4_p.second)* s4 + 8*(c_T3_p.first*c_T4_p.first*alpha0.second.first + s3*c_T4_p.first*c_T3_p.second + s3*c_T3_p.first*c_T4_p.second)* s1*s5*s6)) /(4*pow(3 + alpha0.first,3)*pow(c_T3_p.first,3)*pow(c_T4_p.first,3) -9*s2*s3*c_T3_p.first*c_T4_p.first*s4)))/ (2.*pow(c_T4_p.first,3)*pow(s2,2)*sqrt(-((disc + 2*s1*s5*s6)/(s2*s6)))); double aM_P = (2*(alpha0.first*c_T4_p.first*alpha0.second.second + c_T4_p.second*s2)*(disc + 2*s1*s5*s6) - c_T4_p.first*s2*(4*s3*c_T3_p.first*c_T4_p.first*(c_T3_p.first*c_T4_p.first*alpha0.second.second +s3*c_T4_p.first*c_T3_p.second +s3*c_T3_p.first*c_T4_p.second) + (disc*(- 9*s2*s3*c_T2_p.first*c_T3_p.first*c_T4_p.first* c_T2_p.second + 9*((1 - alpha0.first*(3 + 2*alpha0.first))* c_T3_p.first*c_T4_p.first*alpha0.second.second - s2*s3*c_T4_p.first*c_T3_p.second - s2*s3*c_T3_p.first*c_T4_p.second)* s4 + 8*(c_T3_p.first*c_T4_p.first*alpha0.second.second + s3*c_T4_p.first*c_T3_p.second + s3*c_T3_p.first*c_T4_p.second)* s1*s5*s6)) /(4*pow(3 + alpha0.first,3)*pow(c_T3_p.first,3)*pow(c_T4_p.first,3) -9*s2*s3*c_T3_p.first*c_T4_p.first*s4)))/ (2.*pow(c_T4_p.first,3)*pow(s2,2)*sqrt(-((disc + 2*s1*s5*s6)/(s2*s6)))); double aP_M = (c_T4_p.first*sqrt(-((disc + 2*s1*s5*s6)/ (s2*s6)))* (-2*c_T3_p.second - (2*alpha0.first*c_T3_p.first*alpha0.second.first)/s2 + (c_T3_p.first*(4*s3*c_T3_p.first*c_T4_p.first*(c_T3_p.first*c_T4_p.first*alpha0.second.first + s3*c_T4_p.first*c_T3_p.second + s3*c_T3_p.first*c_T4_p.second) + (disc*(-9*s2*s3*c_T2_p.first*c_T3_p.first*c_T4_p.first* c_T2_p.second + 9*((1 - alpha0.first*(3 + 2*alpha0.first))* c_T3_p.first*c_T4_p.first*alpha0.second.first - s2*s3*c_T4_p.first*c_T3_p.second - s2*s3*c_T3_p.first*c_T4_p.second)* s4 + 8*(c_T3_p.first*c_T4_p.first*alpha0.second.first + s3*c_T4_p.first*c_T3_p.second + s3*c_T3_p.first*c_T4_p.second)* s1*s5*s6))/ (4* pow(3 + alpha0.first,3)* pow(c_T3_p.first,3)* pow(c_T4_p.first,3) - 9*s2*s3*c_T3_p.first*c_T4_p.first*s4)))/ (disc + 2*s1*s5*s6)))/(2.*pow(c_T3_p.first,2)); double aP_P = (c_T4_p.first*sqrt(-((disc + 2*s1*s5*s6)/ (s2*s6)))* (-2*c_T3_p.second - (2*alpha0.first*c_T3_p.first*alpha0.second.second)/s2 + (c_T3_p.first*(4*s3*c_T3_p.first*c_T4_p.first*(c_T3_p.first*c_T4_p.first*alpha0.second.second + s3*c_T4_p.first*c_T3_p.second + s3*c_T3_p.first*c_T4_p.second) + (disc*(-9*s2*s3*c_T2_p.first*c_T3_p.first*c_T4_p.first* c_T2_p.second + 9*((1 - alpha0.first*(3 + 2*alpha0.first))* c_T3_p.first*c_T4_p.first*alpha0.second.second - s2*s3*c_T4_p.first*c_T3_p.second - s2*s3*c_T3_p.first*c_T4_p.second)* s4 + 8*(c_T3_p.first*c_T4_p.first*alpha0.second.second + s3*c_T4_p.first*c_T3_p.second + s3*c_T3_p.first*c_T4_p.second)* s1*s5*s6))/ (4* pow(3 + alpha0.first,3)* pow(c_T3_p.first,3)* pow(c_T4_p.first,3) - 9*s2*s3*c_T3_p.first*c_T4_p.first*s4)))/ (disc + 2*s1*s5*s6)))/(2.*pow(c_T3_p.first,2)); - Scatter2DPtr _h_alphaM = bookScatter2D(1,3,1); + Scatter2DPtr _h_alphaM; + book(_h_alphaM,1,3,1); _h_alphaM->addPoint(0.5, aM, make_pair(0.5,0.5), make_pair(-aM_M , -aM_P ) ); - Scatter2DPtr _h_alphaP = bookScatter2D(1,4,1); + Scatter2DPtr _h_alphaP; + book(_h_alphaP,1,4,1); _h_alphaP->addPoint(0.5, aP, make_pair(0.5,0.5), make_pair(-aP_M , -aP_P ) ); // now for Delta double sDelta = (-2.*(3. + alpha0.first)*c_T3_p.first)/(aM*sqrt(1 - sqr(alpha0.first))); double cDelta = (-3*(3 + alpha0.first)*c_T2_p.first)/(aM*aP*sqrt(1 - sqr(alpha0.first))); double Delta = asin(sDelta); if(cDelta<0.) Delta = M_PI-Delta; double ds_P = (-9*c_T2_p.first*((-1 + alpha0.first)*(1 + alpha0.first)* (3 + alpha0.first)*c_T3_p.first*c_T4_p.first*c_T2_p.second + c_T2_p.first*c_T4_p.first*(c_T3_p.first*(alpha0.second.first + 3*alpha0.first*alpha0.second.first) -(-1 + alpha0.first)*(1 + alpha0.first)*(3 + alpha0.first)*c_T3_p.second) - (-1 + alpha0.first)*(1 + alpha0.first)* (3 + alpha0.first)*c_T2_p.first*c_T3_p.first*c_T4_p.second)*disc)/ (pow(1 - pow(alpha0.first,2),1.5)*pow(c_T4_p.first,3)*pow(-((disc + 2*s1*s5*s6)/ (s2*s6)),1.5)*(-9*s2*s4 + 4*s1*s5*s6)); double ds_M = (-9*c_T2_p.first*((-1 + alpha0.first)*(1 + alpha0.first)* (3 + alpha0.first)*c_T3_p.first*c_T4_p.first*c_T2_p.second + c_T2_p.first*c_T4_p.first*(c_T3_p.first*(alpha0.second.second + 3*alpha0.first*alpha0.second.second) -(-1 + alpha0.first)*(1 + alpha0.first)*(3 + alpha0.first)*c_T3_p.second) - (-1 + alpha0.first)*(1 + alpha0.first)* (3 + alpha0.first)*c_T2_p.first*c_T3_p.first*c_T4_p.second)*disc)/ (pow(1 - pow(alpha0.first,2),1.5)*pow(c_T4_p.first,3)*pow(-((disc + 2*s1*s5*s6)/ (s2*s6)),1.5)*(-9*s2*s4 + 4*s1*s5*s6)); ds_P /= sqrt(1.-sqr(sDelta)); ds_M /= sqrt(1.-sqr(sDelta)); - Scatter2DPtr _h_sin = bookScatter2D(1,2,1); + Scatter2DPtr _h_sin; + book(_h_sin,1,2,1); _h_sin->addPoint(0.5, Delta/M_PI*180., make_pair(0.5,0.5), make_pair( -ds_P/M_PI*180., -ds_M/M_PI*180. ) ); } // alpha 0 pair c_T2_n = calcCoeff(2,_h_T2_n); pair c_T3_n = calcCoeff(3,_h_T3_n); pair c_T4_n = calcCoeff(4,_h_T4_n); s4 = sqr(c_T2_n.first); s5 = sqr(c_T3_n.first); s6 = sqr(c_T4_n.first); disc = s1*s5*s6*(-9.*s2*s4 + 4.*s1*s5*s6); if(disc>=0.) { disc = sqrt(disc); double aM = sqrt(-1./s2/s6*(2.*s1*s5*s6+disc)); double a0 = c_T4_n.first/c_T3_n.first*aM; double a0_M = (c_T4_n.first*sqrt(-((disc + 2*s1*s5*s6)/ (s2*s6)))* (-2*c_T3_n.second - (2*alpha0.first*c_T3_n.first*alpha0.second.first)/s2 + (c_T3_n.first*(4*s3*c_T3_n.first*c_T4_n.first*(c_T3_n.first*c_T4_n.first*alpha0.second.first + s3*c_T4_n.first*c_T3_n.second + s3*c_T3_n.first*c_T4_n.second) + (disc*(-9*s2*s3*c_T2_n.first*c_T3_n.first*c_T4_n.first* c_T2_n.second + 9*((1 - alpha0.first*(3 + 2*alpha0.first))* c_T3_n.first*c_T4_n.first*alpha0.second.first - s2*s3*c_T4_n.first*c_T3_n.second - s2*s3*c_T3_n.first*c_T4_n.second)* s4 + 8*(c_T3_n.first*c_T4_n.first*alpha0.second.first + s3*c_T4_n.first*c_T3_n.second + s3*c_T3_n.first*c_T4_n.second)* s1*s5*s6))/ (4* pow(3 + alpha0.first,3)* pow(c_T3_n.first,3)* pow(c_T4_n.first,3) - 9*s2*s3*c_T3_n.first*c_T4_n.first*s4)))/ (disc + 2*s1*s5*s6)))/(2.*pow(c_T3_n.first,2)); double a0_P = (c_T4_n.first*sqrt(-((disc + 2*s1*s5*s6)/ (s2*s6)))* (-2*c_T3_n.second - (2*alpha0.first*c_T3_n.first*alpha0.second.second)/s2 + (c_T3_n.first*(4*s3*c_T3_n.first*c_T4_n.first*(c_T3_n.first*c_T4_n.first*alpha0.second.second + s3*c_T4_n.first*c_T3_n.second + s3*c_T3_n.first*c_T4_n.second) + (disc*(-9*s2*s3*c_T2_n.first*c_T3_n.first*c_T4_n.first* c_T2_n.second + 9*((1 - alpha0.first*(3 + 2*alpha0.first))* c_T3_n.first*c_T4_n.first*alpha0.second.second - s2*s3*c_T4_n.first*c_T3_n.second - s2*s3*c_T3_n.first*c_T4_n.second)* s4 + 8*(c_T3_n.first*c_T4_n.first*alpha0.second.second + s3*c_T4_n.first*c_T3_n.second + s3*c_T3_n.first*c_T4_n.second)* s1*s5*s6))/ (4* pow(3 + alpha0.first,3)* pow(c_T3_n.first,3)* pow(c_T4_n.first,3) - 9*s2*s3*c_T3_n.first*c_T4_n.first*s4)))/ (disc + 2*s1*s5*s6)))/(2.*pow(c_T3_n.first,2)); - Scatter2DPtr _h_alpha0 = bookScatter2D(1,5,1); + Scatter2DPtr _h_alpha0; + book(_h_alpha0,1,5,1); _h_alpha0->addPoint(0.5, a0, make_pair(0.5,0.5), make_pair(-a0_M , -a0_P ) ); } } //@} /// @name Histograms //@{ Histo1DPtr _h_T1_p,_h_T2_p,_h_T3_p,_h_T4_p,_h_T5_p; Histo1DPtr _h_T1_n,_h_T2_n,_h_T3_n,_h_T4_n,_h_T5_n; Histo1DPtr _h_cThetaL; Histo1DPtr _h_mu_p,_h_mu_n; - double _wsum_p,_wsum_n; + CounterPtr _wsum_p,_wsum_n; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2019_I1691850); } diff --git a/analyses/pluginBES/BESIII_2019_I1723934.cc b/analyses/pluginBES/BESIII_2019_I1723934.cc --- a/analyses/pluginBES/BESIII_2019_I1723934.cc +++ b/analyses/pluginBES/BESIII_2019_I1723934.cc @@ -1,115 +1,116 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2019_I1723934 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2019_I1723934); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); - _nChi0 = bookCounter("TMP/chi0"); + book(_nChi0, "TMP/chi0"); } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - --nRes[child.pdgId()]; + --nRes[child.pid()]; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } const FinalState& ufs = apply(event, "UFS"); bool matched = false; - foreach (const Particle& p, ufs.particles(Cuts::pid==10441)) { + for (const Particle& p : ufs.particles(Cuts::pid==10441)) { if(p.children().empty()) continue; map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); - foreach (const Particle & p2, ufs.particles(Cuts::pid==223)) { + for (const Particle& p2 : ufs.particles(Cuts::pid==223)) { map nRes2 = nRes; int ncount2 = ncount; findChildren(p2,nRes2,ncount2); if(ncount2!=0) continue; matched=true; for(auto const & val : nRes2) { if(val.second!=0) { matched = false; break; } } if(matched) { - _nChi0->fill(event.weight()); + _nChi0->fill(); break; } } if(matched) break; } } /// Normalise histograms etc., after the run void finalize() { double fact = crossSection()/ sumOfWeights() /picobarn; double sigma = _nChi0->val()*fact; double error = _nChi0->err()*fact; Scatter2D temphisto(refData(1, 1, 1)); - Scatter2DPtr mult = bookScatter2D(1, 1, 1); + Scatter2DPtr mult; + book(mult, 1, 1, 1); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/GeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } //@} /// @name Histograms //@{ CounterPtr _nChi0; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2019_I1723934); } diff --git a/analyses/pluginBES/BESIII_2019_I1724880.cc b/analyses/pluginBES/BESIII_2019_I1724880.cc --- a/analyses/pluginBES/BESIII_2019_I1724880.cc +++ b/analyses/pluginBES/BESIII_2019_I1724880.cc @@ -1,71 +1,71 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2019_I1724880 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2019_I1724880); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS"); // Book histograms - _h_n = bookHisto1D(1, 1, 1); + book(_h_n, 1, 1, 1); } void findChildren(const Particle & p,int & nCharged) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - if(PID::isCharged(child.pdgId())) ++nCharged; + if(PID::isCharged(child.pid())) ++nCharged; } else findChildren(child,nCharged); } } /// Perform the per-event analysis void analyze(const Event& event) { - foreach (const Particle& p, apply(event, "UFS").particles(Cuts::pid==441)) { + for (const Particle& p : apply(event, "UFS").particles(Cuts::pid==441)) { int nCharged(0); findChildren(p,nCharged); - _h_n->fill(min(nCharged,8),event.weight()); + _h_n->fill(min(nCharged,8)); } } /// Normalise histograms etc., after the run void finalize() { normalize(_h_n,2.); } //@} /// @name Histograms //@{ Histo1DPtr _h_n; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2019_I1724880); } diff --git a/analyses/pluginBES/BESIII_2019_I1736235.cc b/analyses/pluginBES/BESIII_2019_I1736235.cc --- a/analyses/pluginBES/BESIII_2019_I1736235.cc +++ b/analyses/pluginBES/BESIII_2019_I1736235.cc @@ -1,80 +1,81 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" namespace Rivet { /// @brief Add a short analysis description here class BESIII_2019_I1736235 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2019_I1736235); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(FinalState(), "FS"); // Book histograms - _nproton = bookCounter("TMP/proton"); + book(_nproton, "TMP/proton"); } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); if(fs.particles().size()!=2) vetoEvent; - foreach (const Particle& p, fs.particles()) { - if(abs(p.pdgId())!=PID::PROTON) vetoEvent; + for (const Particle& p : fs.particles()) { + if(abs(p.pid())!=PID::PROTON) vetoEvent; } - _nproton->fill(event.weight()); + _nproton->fill(); } /// Normalise histograms etc., after the run void finalize() { double sigma = _nproton->val(); double error = _nproton->err(); sigma *= crossSection()/ sumOfWeights() /picobarn; error *= crossSection()/ sumOfWeights() /picobarn; Scatter2D temphisto(refData(1, 1, 1)); - Scatter2DPtr mult = bookScatter2D(1, 1, 1); + Scatter2DPtr mult; + book(mult, 1, 1, 1); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/GeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } //@} /// @name Histograms //@{ CounterPtr _nproton; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2019_I1736235); } diff --git a/analyses/pluginBES/BESIII_2019_I1743841.cc b/analyses/pluginBES/BESIII_2019_I1743841.cc --- a/analyses/pluginBES/BESIII_2019_I1743841.cc +++ b/analyses/pluginBES/BESIII_2019_I1743841.cc @@ -1,129 +1,130 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief 2K+2K- and K+K- phi cross section class BESIII_2019_I1743841 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(BESIII_2019_I1743841); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); // book counters - _c2Kp2Km = bookCounter("TMP/2Kp2Km" ); - _cKpKmPhi = bookCounter("TMP/KpKmPhi"); + book(_c2Kp2Km , "TMP/2Kp2Km" ); + book(_cKpKmPhi, "TMP/KpKmPhi"); } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - nRes[child.pdgId()]-=1; + nRes[child.pid()]-=1; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } if(ntotal==4 && nCount[321]==2 && nCount[-321]==2) - _c2Kp2Km->fill(event.weight()); + _c2Kp2Km->fill(); const FinalState& ufs = apply(event, "UFS"); - foreach (const Particle& p, ufs.particles(Cuts::pid==333)) { + for (const Particle& p : ufs.particles(Cuts::pid==333)) { if(p.children().empty()) continue; map nRes=nCount; int ncount = ntotal; findChildren(p,nRes,ncount); // phi K+K- if(ncount==2) { bool matched = true; for(auto const & val : nRes) { if(abs(val.first)==321) { if(val.second!=1) { matched = false; break; } } else if(val.second!=0) { matched = false; break; } } if(matched) - _cKpKmPhi->fill(event.weight()); + _cKpKmPhi->fill(); } } } /// Normalise histograms etc., after the run void finalize() { for(unsigned int ix=1;ix<3;++ix) { double sigma = 0., error = 0.; if(ix==1) { sigma = _c2Kp2Km->val(); error = _c2Kp2Km->err(); } else { sigma = _cKpKmPhi->val(); error = _cKpKmPhi->err(); } sigma *= crossSection()/ sumOfWeights() /picobarn; error *= crossSection()/ sumOfWeights() /picobarn; Scatter2D temphisto(refData(ix, 1, 1)); - Scatter2DPtr mult = bookScatter2D(ix, 1, 1); + Scatter2DPtr mult; + book(mult,ix, 1, 1); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/GeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } } //@} /// @name Histograms //@{ CounterPtr _c2Kp2Km,_cKpKmPhi; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(BESIII_2019_I1743841); } diff --git a/analyses/pluginCESR/CLEO_2000_I526554.cc b/analyses/pluginCESR/CLEO_2000_I526554.cc --- a/analyses/pluginCESR/CLEO_2000_I526554.cc +++ b/analyses/pluginCESR/CLEO_2000_I526554.cc @@ -1,116 +1,114 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/Beam.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class CLEO_2000_I526554 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(CLEO_2000_I526554); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { declare(Beam(), "Beams"); declare(UnstableParticles(), "UFS"); // Book histograms - _h_Ds_star1 = bookHisto1D(1, 1, 1); - _h_Ds = bookHisto1D(2, 1, 1); - _h_Ds_star2 = bookHisto1D(3, 1, 1); - _h_Ds_primary = bookHisto1D(4, 1, 1); + book(_h_Ds_star1 , 1, 1, 1); + book(_h_Ds , 2, 1, 1); + book(_h_Ds_star2 , 3, 1, 1); + book(_h_Ds_primary, 4, 1, 1); } /// Perform the per-event analysis void analyze(const Event& event) { - const double weight = event.weight(); - // Loop through unstable FS particles and look for charmed mesons const UnstableParticles& ufs = apply(event, "UFS"); const Beam beamproj = apply(event, "Beams"); const ParticlePair& beams = beamproj.beams(); const FourMomentum mom_tot = beams.first.momentum() + beams.second.momentum(); LorentzTransform cms_boost; if (mom_tot.p3().mod() > 1*MeV) cms_boost = LorentzTransform::mkFrameTransformFromBeta(mom_tot.betaVec()); const double s = sqr(beamproj.sqrtS()); for (const Particle& p : ufs.particles(Cuts::abspid==431 or Cuts::abspid==433)) { // 3-momentum in CMS frame const double mom = cms_boost.transform(p.momentum()).vector3().mod(); const int pdgid = p.abspid(); double mH2(0.),xp(0.); bool primary = true; switch (pdgid) { case 431: // MSG_DEBUG("D_s found"); mH2 = sqr(1.96834); xp = mom/sqrt(s/4.0 - mH2); - _h_Ds->fill(xp,weight); + _h_Ds->fill(xp); for(const Particle & mother : p.parents()) { - if(PID::isCharmMeson(mother.pdgId())) { + if(PID::isCharmMeson(mother.pid())) { primary = false; break; } } if(primary) - _h_Ds_primary->fill(xp,weight); + _h_Ds_primary->fill(xp); break; case 433: MSG_DEBUG("D_s* found"); mH2 = sqr(2.1122); xp = mom/sqrt(s/4.0 - mH2); - _h_Ds_star1->fill(xp,weight); - _h_Ds_star2->fill(xp,weight); + _h_Ds_star1->fill(xp); + _h_Ds_star2->fill(xp); break; default: break; } } } /// Normalise histograms etc., after the run void finalize() { // BR(D_s->phi pi) x BR(phi->K+K-) double br2 = 0.0227; // x D_s* -> gamma d_s double br1 = br2*.935; // cross section factor double fact = crossSection()/picobarn/sumOfWeights(); // normalize the cross sections (bin width) scale(_h_Ds_star1 , br1*fact*0.03 ); scale(_h_Ds , br2*fact*0.03 ); scale(_h_Ds_star2 , br2*fact*0.03 ); scale(_h_Ds_primary, br2*fact*0.03 ); } //@} /// @name Histograms //@{ Histo1DPtr _h_Ds,_h_Ds_primary,_h_Ds_star1,_h_Ds_star2; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(CLEO_2000_I526554); } diff --git a/analyses/pluginCESR/CLEO_2000_I537236.cc b/analyses/pluginCESR/CLEO_2000_I537236.cc --- a/analyses/pluginCESR/CLEO_2000_I537236.cc +++ b/analyses/pluginCESR/CLEO_2000_I537236.cc @@ -1,113 +1,114 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Xi_c0 -> Xi-pi+ asymmetry class CLEO_2000_I537236 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(CLEO_2000_I537236); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS" ); // Book histograms - _h_ctheta = bookHisto1D("ctheta" , 20,-1,1); + book(_h_ctheta, "ctheta" , 20,-1,1); } /// Perform the per-event analysis void analyze(const Event& event) { // loop over Omega baryons - foreach(const Particle& Xic, apply(event, "UFS").particles(Cuts::abspid==4132)) { - int sign = Xic.pdgId()/4132; + for( const Particle& Xic : apply(event, "UFS").particles(Cuts::abspid==4132)) { + int sign = Xic.pid()/4132; if(Xic.children().size()!=2) continue; Particle baryon1,meson1; - if(Xic.children()[0].pdgId()==sign*3312 && - Xic.children()[1].pdgId()==sign*211) { + if(Xic.children()[0].pid()==sign*3312 && + Xic.children()[1].pid()==sign*211) { baryon1 = Xic.children()[0]; meson1 = Xic.children()[1]; } - else if(Xic.children()[1].pdgId()==sign*3312 && - Xic.children()[0].pdgId()==sign*211) { + else if(Xic.children()[1].pid()==sign*3312 && + Xic.children()[0].pid()==sign*211) { baryon1 = Xic.children()[1]; meson1 = Xic.children()[0]; } else continue; Particle baryon2,meson2; - if(baryon1.children()[0].pdgId()== sign*3122 && - baryon1.children()[1].pdgId()==-sign*211) { + if(baryon1.children()[0].pid()== sign*3122 && + baryon1.children()[1].pid()==-sign*211) { baryon2 = baryon1.children()[0]; meson2 = baryon1.children()[1]; } - else if(baryon1.children()[1].pdgId()== sign*3122 && - baryon1.children()[0].pdgId()==-sign*211) { + else if(baryon1.children()[1].pid()== sign*3122 && + baryon1.children()[0].pid()==-sign*211) { baryon2 = baryon1.children()[1]; meson2 = baryon1.children()[0]; } else continue; // first boost to the Xic rest frame LorentzTransform boost1 = LorentzTransform::mkFrameTransformFromBeta(Xic.momentum().betaVec()); FourMomentum pbaryon1 = boost1.transform(baryon1.momentum()); FourMomentum pbaryon2 = boost1.transform(baryon2.momentum()); // to lambda rest frame LorentzTransform boost2 = LorentzTransform::mkFrameTransformFromBeta(pbaryon1.betaVec()); Vector3 axis = pbaryon1.p3().unit(); FourMomentum pp = boost2.transform(pbaryon2); // calculate angle double cTheta = pp.p3().unit().dot(axis); _h_ctheta->fill(cTheta,1.); } } pair calcAlpha(Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,0.); double sum1(0.),sum2(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double ai = 0.5*(bin.xMax()-bin.xMin()); double bi = 0.5*ai*(bin.xMax()+bin.xMin()); double Ei = bin.areaErr(); sum1 += sqr(bi/Ei); sum2 += bi/sqr(Ei)*(Oi-ai); } return make_pair(sum2/sum1,sqrt(1./sum1)); } /// Normalise histograms etc., after the run void finalize() { normalize(_h_ctheta); - Scatter2DPtr _h_alpha = bookScatter2D(1,1,1); + Scatter2DPtr _h_alpha; + book(_h_ctheta,1,1,1); pair alpha = calcAlpha(_h_ctheta); _h_alpha->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second,alpha.second) ); } //@} /// @name Histograms //@{ Histo1DPtr _h_ctheta; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(CLEO_2000_I537236); } diff --git a/analyses/pluginCESR/CLEO_2001_I552541.cc b/analyses/pluginCESR/CLEO_2001_I552541.cc --- a/analyses/pluginCESR/CLEO_2001_I552541.cc +++ b/analyses/pluginCESR/CLEO_2001_I552541.cc @@ -1,181 +1,181 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class CLEO_2001_I552541 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(CLEO_2001_I552541); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS"); - _d_Dbar0[0] = bookCounter("/TMP/d_D0_low" ); - _d_Dbar0[1] = bookCounter("/TMP/d_D0_high"); - _d_Dm[0] = bookCounter("/TMP/d_Dm_low" ); - _d_Dm[1] = bookCounter("/TMP/d_Dm_high"); - _d_Lam[0] = bookCounter("/TMP/d_La_low" ); - _d_Lam[1] = bookCounter("/TMP/d_La_high"); + book(_d_Dbar0[0], "/TMP/d_D0_low" ); + book(_d_Dbar0[1], "/TMP/d_D0_high"); + book(_d_Dm[0] , "/TMP/d_Dm_low" ); + book(_d_Dm[1] , "/TMP/d_Dm_high"); + book(_d_Lam[0] , "/TMP/d_La_low" ); + book(_d_Lam[1] , "/TMP/d_La_high"); - _n_Dbar0[0][0] = bookCounter("/TMP/d_D0_low_low" ); - _n_Dbar0[0][1] = bookCounter("/TMP/d_D0_low_high" ); - _n_Dbar0[1][0] = bookCounter("/TMP/d_D0_high_low" ); - _n_Dbar0[1][1] = bookCounter("/TMP/d_D0_high_high"); - _n_Dm[0][0] = bookCounter("/TMP/d_Dm_low_low" ); - _n_Dm[0][1] = bookCounter("/TMP/d_Dm_low_high" ); - _n_Dm[1][0] = bookCounter("/TMP/d_Dm_high_low" ); - _n_Dm[1][1] = bookCounter("/TMP/d_Dm_high_high"); - _n_Lam[0][0] = bookCounter("/TMP/d_La_low_low" ); - _n_Lam[0][1] = bookCounter("/TMP/d_La_low_high" ); - _n_Lam[1][0] = bookCounter("/TMP/d_La_high_low" ); - _n_Lam[1][1] = bookCounter("/TMP/d_La_high_high"); + book(_n_Dbar0[0][0], "/TMP/d_D0_low_low" ); + book(_n_Dbar0[0][1], "/TMP/d_D0_low_high" ); + book(_n_Dbar0[1][0], "/TMP/d_D0_high_low" ); + book(_n_Dbar0[1][1], "/TMP/d_D0_high_high"); + book(_n_Dm[0][0] , "/TMP/d_Dm_low_low" ); + book(_n_Dm[0][1] , "/TMP/d_Dm_low_high" ); + book(_n_Dm[1][0] , "/TMP/d_Dm_high_low" ); + book(_n_Dm[1][1] , "/TMP/d_Dm_high_high"); + book(_n_Lam[0][0] , "/TMP/d_La_low_low" ); + book(_n_Lam[0][1] , "/TMP/d_La_low_high" ); + book(_n_Lam[1][0] , "/TMP/d_La_high_low" ); + book(_n_Lam[1][1] , "/TMP/d_La_high_high"); } /// Perform the per-event analysis void analyze(const Event& event) { - double weight=event.weight(); const UnstableParticles& ufs = apply(event, "UFS"); for(const Particle & p : ufs.particles(Cuts::pid==-4122 or Cuts::pid==-411 or Cuts::pid==-421)) { - long id1 = p.pdgId(); + long id1 = p.pid(); double mom1 = p.p3().mod(); if(mom1<2.3*GeV || mom1>5.*GeV) continue; bool high1 = mom1>3.3*GeV; if(id1==-4122) { - _d_Lam[high1]->fill(weight); + _d_Lam[high1]->fill(); } else if(id1==-411) { - _d_Dm[high1]->fill(weight); + _d_Dm[high1]->fill(); } else if(id1==-421) { - _d_Dbar0[high1]->fill(weight); + _d_Dbar0[high1]->fill(); } for(const Particle & p2 : ufs.particles(Cuts::pid==4122)) { if(p.p3().angle(p2.p3())<0.5*M_PI) continue; double mom2 = p2.p3().mod(); if(mom2<2.3*GeV || mom2>5.*GeV) continue; bool high2 = mom2>3.3*GeV; if(id1==-4122) { - _n_Lam[high1][high2]->fill(weight); + _n_Lam[high1][high2]->fill(); } else if(id1==-411) { - _n_Dm[high1][high2]->fill(weight); + _n_Dm[high1][high2]->fill(); } else if(id1==-421) { - _n_Dbar0[high1][high2]->fill(weight); + _n_Dbar0[high1][high2]->fill(); } } } } /// Normalise histograms etc., after the run void finalize() { Scatter1D R_D0_low_low = *_n_Dbar0[0][0]/ *_d_Dbar0[0]; Scatter1D R_D0_low_high = *_n_Dbar0[0][1]/ *_d_Dbar0[0]; Scatter1D R_D0_high_low = *_n_Dbar0[1][0]/ *_d_Dbar0[1]; Scatter1D R_D0_high_high = *_n_Dbar0[1][1]/ *_d_Dbar0[1]; Scatter1D R_Dm_low_low = *_n_Dm[0][0] / *_d_Dm[0]; Scatter1D R_Dm_low_high = *_n_Dm[0][1] / *_d_Dm[0]; Scatter1D R_Dm_high_low = *_n_Dm[1][0] / *_d_Dm[1]; Scatter1D R_Dm_high_high = *_n_Dm[1][1] / *_d_Dm[1]; Scatter1D R_Lam_low_low = *_n_Lam[0][0] / *_d_Lam[0]; Scatter1D R_Lam_low_high = *_n_Lam[0][1] / *_d_Lam[0]; Scatter1D R_Lam_high_low = *_n_Lam[1][0] / *_d_Lam[1]; Scatter1D R_Lam_high_high = *_n_Lam[1][1] / *_d_Lam[1]; for(unsigned int ix=3;ix<5;++ix) { for(unsigned int iy=1;iy<5;++iy) { double num(0.),den(0.),num_err(0.),den_err(0.); if(ix==3) { if(iy==1) { den = R_D0_low_low .points()[0].x(); den_err = R_D0_low_low .points()[0].xErrAvg(); } else if(iy==2) { den = R_D0_high_low .points()[0].x(); den_err = R_D0_high_low .points()[0].xErrAvg(); } else if(iy==3) { den = R_D0_low_low .points()[0].x(); den_err = R_D0_low_low .points()[0].xErrAvg(); } else if(iy==4) { den = R_D0_high_high.points()[0].x(); den_err = R_D0_high_high.points()[0].xErrAvg(); } } else if(ix==4) { if(iy==1) { den = R_Dm_low_low .points()[0].x(); den_err = R_Dm_low_low .points()[0].xErrAvg(); } else if(iy==2) { den = R_Dm_high_low .points()[0].x(); den_err = R_Dm_high_low .points()[0].xErrAvg(); } else if(iy==3) { den = R_Dm_low_low .points()[0].x(); den_err = R_Dm_low_low .points()[0].xErrAvg(); } else if(iy==4) { den = R_Dm_high_high.points()[0].x(); den_err = R_Dm_high_high.points()[0].xErrAvg(); } } if(iy==1) { num = R_Lam_low_low .points()[0].x(); num_err = R_Lam_low_low .points()[0].xErrAvg(); } else if(iy==2) { num = R_Lam_high_low .points()[0].x(); num_err = R_Lam_high_low .points()[0].xErrAvg(); } else if(iy==3) { num = R_Lam_low_low .points()[0].x(); num_err = R_Lam_low_low .points()[0].xErrAvg(); } else if(iy==4) { num = R_Lam_high_high.points()[0].x(); num_err = R_Lam_high_high.points()[0].xErrAvg(); } double val = num/den; double err = val>=0. ? val*sqrt(sqr(num_err/num)+sqr(den_err/den)) : 0.; - Scatter2DPtr ratio = bookScatter2D(ix, 1, iy); + Scatter2DPtr ratio; + book(ratio,ix, 1, iy); Scatter2D temphisto(refData(ix, 1, iy)); const double x = temphisto.point(0).x(); pair ex = temphisto.point(0).xErrs(); ratio->addPoint(x, val, ex, make_pair(err,err)); } } } //@} /// @name Histograms //@{ CounterPtr _d_Dbar0[2],_d_Dm[2],_d_Lam[2]; CounterPtr _n_Dbar0[2][2],_n_Dm[2][2],_n_Lam[2][2]; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(CLEO_2001_I552541); } diff --git a/analyses/pluginFrascati/FENICE_1996_I426675.cc b/analyses/pluginFrascati/FENICE_1996_I426675.cc --- a/analyses/pluginFrascati/FENICE_1996_I426675.cc +++ b/analyses/pluginFrascati/FENICE_1996_I426675.cc @@ -1,86 +1,86 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" namespace Rivet { /// @brief Add a short analysis description here class FENICE_1996_I426675 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(FENICE_1996_I426675); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(FinalState(), "FS"); // Book histograms - _c_hadrons = bookCounter("/TMP/sigma_hadrons"); - + book(_c_hadrons, "/TMP/sigma_hadrons"); } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } // mu+mu- + photons if(nCount[-13]==1 and nCount[13]==1 && ntotal==2+nCount[22]) vetoEvent; - _c_hadrons->fill(event.weight()); + _c_hadrons->fill(); } /// Normalise histograms etc., after the run void finalize() { double fact = crossSection()/ sumOfWeights() /nanobarn; double sig_h = _c_hadrons->val()*fact; double err_h = _c_hadrons->err()*fact; Scatter2D temphisto(refData(1, 1, 1)); - Scatter2DPtr hadrons = bookScatter2D(1, 1, 1); + Scatter2DPtr hadrons; + book(hadrons, 1, 1, 1); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqr(sqrtS()/GeV), x-ex2.first, x+ex2.second)) { hadrons->addPoint(x, sig_h, ex, make_pair(err_h,err_h)); } else { hadrons->addPoint(x, 0., ex, make_pair(0.,.0)); } } } //@} /// @name Histograms //@{ CounterPtr _c_hadrons; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(FENICE_1996_I426675); } diff --git a/analyses/pluginFrascati/KLOE2_2014_I1317236.cc b/analyses/pluginFrascati/KLOE2_2014_I1317236.cc --- a/analyses/pluginFrascati/KLOE2_2014_I1317236.cc +++ b/analyses/pluginFrascati/KLOE2_2014_I1317236.cc @@ -1,108 +1,108 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief form factor for phi-> eta class KLOE2_2014_I1317236 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(KLOE2_2014_I1317236); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS"); // Book histograms - _h_m = bookHisto1D(1, 1, 1); - _weight=0.; + book(_h_m, 1, 1, 1); + book(_weight,"TMP/weight"); } void findDecayProducts(const Particle & mother, unsigned int & nstable, unsigned int & neta, unsigned int & nep, unsigned int & nem, unsigned int & ngamma, FourMomentum & ptot) { for(const Particle & p : mother.children()) { - int id = p.pdgId(); + int id = p.pid(); if (id == PID::EMINUS ) { ++nem; ++nstable; ptot += p.momentum(); } else if (id == PID::EPLUS) { ++nep; ++nstable; ptot += p.momentum(); } else if (id == PID::ETA) { ++neta; ++nstable; } else if (id == PID::GAMMA) { ++ngamma; ++nstable; } else if ( !p.children().empty() ) { findDecayProducts(p, nstable, neta,nep,nem,ngamma,ptot); } else ++nstable; } } /// Perform the per-event analysis void analyze(const Event& event) { static double me = 0.5109989461*MeV; static double mphi = 1019.461*MeV; static double meta = 547.862*MeV; // Loop over phi mesons - foreach(const Particle& p, apply(event, "UFS").particles(Cuts::pid==333)) { + for(const Particle& p : apply(event, "UFS").particles(Cuts::pid==333)) { unsigned nstable(0),neta(0),nep(0),nem(0),ngamma(0); FourMomentum ptot; findDecayProducts(p,nstable,neta,nep,nem,ngamma,ptot); if(nstable==3 && nem==1 && nem==1 && neta==1) { double q = ptot.mass(); double beta = sqrt(1.-4*sqr(me/q)); double p = sqrt(sqr(1.+sqr(q)/(sqr(mphi)-sqr(meta)))-4.*sqr(mphi*q/(sqr(mphi)-sqr(meta)))); double fact = beta*MeV/q*(1.+2.*sqr(me/q))*pow(p,3); - _h_m->fill(q/MeV,event.weight()/fact); + _h_m->fill(q/MeV,1./fact); } else if(nstable==2 && ngamma ==1 && neta==1) { - _weight += event.weight(); + _weight->fill(); } } } /// Normalise histograms etc., after the run void finalize() { static double alpha= 7.2973525664e-3; - scale(_h_m, 1.5*M_PI/_weight/alpha); + scale(_h_m, 1.5*M_PI/alpha/ *_weight); } //@} /// @name Histograms //@{ Histo1DPtr _h_m; - double _weight; + CounterPtr _weight; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(KLOE2_2014_I1317236); } diff --git a/analyses/pluginFrascati/KLOE2_2016_I1416825.cc b/analyses/pluginFrascati/KLOE2_2016_I1416825.cc --- a/analyses/pluginFrascati/KLOE2_2016_I1416825.cc +++ b/analyses/pluginFrascati/KLOE2_2016_I1416825.cc @@ -1,109 +1,109 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief form factor for phi-> pi0 class KLOE2_2016_I1416825 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(KLOE2_2016_I1416825); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS"); // Book histograms - _h_m = bookHisto1D(1, 1, 1); - _weight=0.; + book(_h_m, 1, 1, 1); + book(_weight,"TMP/weight"); } void findDecayProducts(const Particle & mother, unsigned int & nstable, unsigned int & npi, unsigned int & nep, unsigned int & nem, unsigned int & ngamma, FourMomentum & ptot) { for(const Particle & p : mother.children()) { - int id = p.pdgId(); + int id = p.pid(); if (id == PID::EMINUS ) { ++nem; ++nstable; ptot += p.momentum(); } else if (id == PID::EPLUS) { ++nep; ++nstable; ptot += p.momentum(); } else if (id == PID::PI0) { ++npi; ++nstable; } else if (id == PID::GAMMA) { ++ngamma; ++nstable; } else if ( !p.children().empty() ) { findDecayProducts(p, nstable, npi,nep,nem,ngamma,ptot); } else ++nstable; } } /// Perform the per-event analysis void analyze(const Event& event) { static double me = 0.5109989461*MeV; static double mphi = 1019.461*MeV; static double mpi = 134.9770*MeV; // Loop over phi mesons - foreach(const Particle& p, apply(event, "UFS").particles(Cuts::pid==333)) { + for(const Particle& p : apply(event, "UFS").particles(Cuts::pid==333)) { unsigned nstable(0),npi(0),nep(0),nem(0),ngamma(0); FourMomentum ptot; findDecayProducts(p,nstable,npi,nep,nem,ngamma,ptot); if(nstable==3 && nem==1 && nem==1 && npi==1) { double q = ptot.mass(); double beta = sqrt(1.-4*sqr(me/q)); double p = sqrt(sqr(1.+sqr(q)/(sqr(mphi)-sqr(mpi)))-4.*sqr(mphi*q/(sqr(mphi)-sqr(mpi)))); double fact = beta*MeV/q*(1.+2.*sqr(me/q))*pow(p,3); - _h_m->fill(q/MeV,event.weight()/fact); + _h_m->fill(q/MeV,1./fact); } else if(nstable==2 && ngamma ==1 && npi==1) { - _weight += event.weight(); + _weight->fill(); } } } /// Normalise histograms etc., after the run void finalize() { static double alpha= 7.2973525664e-3; - scale(_h_m, 1.5*M_PI/_weight/alpha); + scale(_h_m, 1.5*M_PI/alpha/ *_weight); } //@} /// @name Histograms //@{ Histo1DPtr _h_m; - double _weight; + CounterPtr _weight; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(KLOE2_2016_I1416825); } diff --git a/analyses/pluginFrascati/KLOE2_2016_I1416990.cc b/analyses/pluginFrascati/KLOE2_2016_I1416990.cc --- a/analyses/pluginFrascati/KLOE2_2016_I1416990.cc +++ b/analyses/pluginFrascati/KLOE2_2016_I1416990.cc @@ -1,125 +1,115 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief eta -> 3 pi analysis class KLOE2_2016_I1416990 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(KLOE2_2016_I1416990); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS"); // Book histograms - _dalitz = bookHisto2D(1,1,1); - _h_dalitz.addHistogram( -0.9, -0.8, bookHisto1D(2, 1, 1 )); - _h_dalitz.addHistogram( -0.8, -0.7, bookHisto1D(2, 1, 2 )); - _h_dalitz.addHistogram( -0.7, -0.6, bookHisto1D(2, 1, 3 )); - _h_dalitz.addHistogram( -0.6, -0.5, bookHisto1D(2, 1, 4 )); - _h_dalitz.addHistogram( -0.5, -0.4, bookHisto1D(2, 1, 5 )); - _h_dalitz.addHistogram( -0.4, -0.3, bookHisto1D(2, 1, 6 )); - _h_dalitz.addHistogram( -0.3, -0.2, bookHisto1D(2, 1, 7 )); - _h_dalitz.addHistogram( -0.2, -0.1, bookHisto1D(2, 1, 8 )); - _h_dalitz.addHistogram( -0.1, 0.0, bookHisto1D(2, 1, 9 )); - _h_dalitz.addHistogram( 0.0, 0.1, bookHisto1D(2, 1, 10)); - _h_dalitz.addHistogram( 0.1, 0.2, bookHisto1D(2, 1, 11)); - _h_dalitz.addHistogram( 0.2, 0.3, bookHisto1D(2, 1, 12)); - _h_dalitz.addHistogram( 0.3, 0.4, bookHisto1D(2, 1, 13)); - _h_dalitz.addHistogram( 0.4, 0.5, bookHisto1D(2, 1, 14)); - _h_dalitz.addHistogram( 0.5, 0.6, bookHisto1D(2, 1, 15)); - _h_dalitz.addHistogram( 0.6, 0.7, bookHisto1D(2, 1, 16)); - _h_dalitz.addHistogram( 0.7, 0.8, bookHisto1D(2, 1, 17)); - _norm=0.; + book(_dalitz, 1,1,1); + double step(0.1), x(-0.9); + for(unsigned int ix=1;ix<18;++ix) { + Histo1DPtr temp; + book(temp,2,1,ix); + _h_dalitz.add( x, x+step, temp); + x += step; + } + book(_norm,"TMP/norm"); } void findDecayProducts(const Particle & mother, unsigned int & nstable, Particles &pi0, Particles &pip, Particles &pim) { for(const Particle & p : mother.children()) { - int id = p.pdgId(); + int id = p.pid(); if (id == PID::PIMINUS ) { pim.push_back(p); ++nstable; } else if (id == PID::PIPLUS) { pip.push_back(p); ++nstable; } else if (id == PID::PI0) { pi0.push_back(p); ++nstable; } else if ( !p.children().empty() ) { findDecayProducts(p, nstable, pi0,pip,pim); } else ++nstable; } } /// Perform the per-event analysis void analyze(const Event& event) { // Loop over eta mesons - foreach(const Particle& p, apply(event, "UFS").particles()) { - if(p.pdgId()!=221) continue; + for (const Particle& p : apply(event, "UFS").particles()) { + if(p.pid()!=221) continue; Particles pi0, pip, pim; unsigned nstable(0); findDecayProducts(p,nstable,pi0,pip,pim); if(nstable==3 && pi0.size()==1 && pip.size()==1 && pim.size()==1) { // massesx double meta = p.mass(); double mpip = pip[0].mass(); double mpim = pim[0].mass(); double mpi0 = pi0[0].mass(); // kinetic energies double Q = meta-mpip-mpim-mpi0; double Ep = 0.5/meta*(sqr(meta)+sqr(mpip)-(p.momentum()-pip[0].momentum()).mass2())-mpip; double Em = 0.5/meta*(sqr(meta)+sqr(mpim)-(p.momentum()-pim[0].momentum()).mass2())-mpim; double E0 = 0.5/meta*(sqr(meta)+sqr(mpi0)-(p.momentum()-pi0[0].momentum()).mass2())-mpi0; // X and Y variables double X = sqrt(3.)/Q*(Ep-Em); double Y = 3.*E0/Q-1.; - _dalitz ->fill(X,Y,event.weight()); - _h_dalitz.fill(Y,X,event.weight()); - if(fabs(X)<0.03225806451612903 && Y>0. && Y<0.1) _norm+=event.weight(); + _dalitz ->fill(X,Y); + _h_dalitz.fill(Y,X); + if(fabs(X)<0.03225806451612903 && Y>0. && Y<0.1) _norm->fill(); } } } /// Normalise histograms etc., after the run void finalize() { - scale(_dalitz,0.06451612903225806*0.1/_norm); - for (Histo1DPtr histo : _h_dalitz.getHistograms()) scale(histo,0.06451612903225806/_norm); + scale(_dalitz,0.06451612903225806*0.1/ *_norm); + for (Histo1DPtr histo : _h_dalitz.histos()) scale(histo,0.06451612903225806/ *_norm); } //@} /// @name Histograms //@{ Histo2DPtr _dalitz; - BinnedHistogram _h_dalitz; - double _norm; + BinnedHistogram _h_dalitz; + CounterPtr _norm; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(KLOE2_2016_I1416990); } diff --git a/analyses/pluginFrascati/KLOE2_2017_I1634981.cc b/analyses/pluginFrascati/KLOE2_2017_I1634981.cc --- a/analyses/pluginFrascati/KLOE2_2017_I1634981.cc +++ b/analyses/pluginFrascati/KLOE2_2017_I1634981.cc @@ -1,78 +1,79 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/FastJets.hh" namespace Rivet { /// @brief Add a short analysis description here class KLOE2_2017_I1634981 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(KLOE2_2017_I1634981); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(FinalState(), "FS"); // Book histograms - _npion = bookCounter("TMP/pion"); + book(_npion, "TMP/pion"); } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); if(fs.particles().size()!=2) vetoEvent; - foreach (const Particle& p, fs.particles()) { - if(abs(p.pdgId())!=PID::PIPLUS) vetoEvent; + for (const Particle& p : fs.particles()) { + if(abs(p.pid())!=PID::PIPLUS) vetoEvent; } - _npion->fill(event.weight()); + _npion->fill(); } /// Normalise histograms etc., after the run void finalize() { double sigma = _npion->val(); double error = _npion->err(); sigma *= crossSection()/ sumOfWeights() /nanobarn; error *= crossSection()/ sumOfWeights() /nanobarn; Scatter2D temphisto(refData(1, 1, 1)); - Scatter2DPtr mult = bookScatter2D(1, 1, 1); + Scatter2DPtr mult; + book(mult, 1, 1, 1); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqr(sqrtS()/GeV), x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } //@} /// @name Histograms //@{ CounterPtr _npion; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(KLOE2_2017_I1634981); } diff --git a/analyses/pluginFrascati/KLOE_2002_I585183.cc b/analyses/pluginFrascati/KLOE_2002_I585183.cc --- a/analyses/pluginFrascati/KLOE_2002_I585183.cc +++ b/analyses/pluginFrascati/KLOE_2002_I585183.cc @@ -1,93 +1,93 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief phi -> pi0 pi0 gamma class KLOE_2002_I585183 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(KLOE_2002_I585183); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { declare(UnstableParticles(), "UFS"); - _h_pipi = bookHisto1D( 1, 1, 1); - _nPhi=0.; + book(_h_pipi, 1, 1, 1); + book(_nPhi,"TMP/nPhi"); } void findDecayProducts(const Particle & mother, unsigned int & nstable, unsigned int & npi, unsigned int & ngamma, FourMomentum & ptot) { for(const Particle & p : mother.children()) { - int id = p.pdgId(); + int id = p.pid(); if (id == PID::PI0) { ++npi; ++nstable; ptot += p.momentum(); } else if (id == PID::GAMMA) { ++ngamma; ++nstable; } else if (id == PID::PIPLUS || id == PID::PIMINUS) { ++nstable; } else if ( !p.children().empty() ) { findDecayProducts(p, nstable, npi, ngamma, ptot); } else ++nstable; } } /// Perform the per-event analysis void analyze(const Event& event) { // Loop over phis for(const Particle& phi : apply(event, "UFS").particles(Cuts::abspid==PID::PHI)) { - _nPhi+=event.weight(); + _nPhi->fill(); unsigned int nstable(0),npi(0),ngamma(0); FourMomentum p_tot(0,0,0,0); findDecayProducts(phi, nstable, npi, ngamma, p_tot); if(nstable!=3) continue; if(npi==2 && ngamma==1 ) { - _h_pipi->fill(p_tot.mass()/MeV, event.weight()); + _h_pipi->fill(p_tot.mass()/MeV); } } } /// Normalise histograms etc., after the run void finalize() { - scale( _h_pipi, 1./_nPhi); + scale( _h_pipi, 1./ *_nPhi); } //@} /// @name Histograms //@{ Histo1DPtr _h_pipi; - double _nPhi; + CounterPtr _nPhi; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(KLOE_2002_I585183); } diff --git a/analyses/pluginMisc/A2_2017_I1486671.cc b/analyses/pluginMisc/A2_2017_I1486671.cc --- a/analyses/pluginMisc/A2_2017_I1486671.cc +++ b/analyses/pluginMisc/A2_2017_I1486671.cc @@ -1,126 +1,126 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief form factors for omega->pi and eta->gamma class A2_2017_I1486671 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(A2_2017_I1486671); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS"); // Book histograms - _h_eta = bookHisto1D(1, 1, 1); - _h_omega = bookHisto1D(2, 1, 1); - _weight_eta =0.; - _weight_omega=0.; + book(_h_eta , 1, 1, 1); + book(_h_omega, 2, 1, 1); + book(_weight_eta ,"TMP/weight_eta"); + book(_weight_omega,"TMP/weight_omega"); } void findDecayProducts(const Particle & mother, unsigned int & nstable, unsigned int & npi, unsigned int & nep, unsigned int & nem, unsigned int & ngamma, FourMomentum & ptot) { for(const Particle & p : mother.children()) { - int id = p.pdgId(); + int id = p.pid(); if (id == PID::EMINUS ) { ++nem; ++nstable; ptot += p.momentum(); } else if (id == PID::EPLUS) { ++nep; ++nstable; ptot += p.momentum(); } else if (id == PID::PI0) { ++npi; ++nstable; } else if (id == PID::GAMMA) { ++ngamma; ++nstable; } else if ( !p.children().empty() ) { findDecayProducts(p, nstable, npi,nep,nem,ngamma,ptot); } else ++nstable; } } /// Perform the per-event analysis void analyze(const Event& event) { static double me = 0.5109989461*MeV; static double momega = 1019.461*MeV; static double meta = 547.862 *MeV; static double mpi = 134.9770*MeV; // Loop over eta and omega mesons - foreach(const Particle& p, apply(event, "UFS").particles(Cuts::pid==333 or Cuts::pid==221)) { + for( const Particle& p : apply(event, "UFS").particles(Cuts::pid==333 or Cuts::pid==221)) { unsigned nstable(0),npi(0),nep(0),nem(0),ngamma(0); FourMomentum ptot; findDecayProducts(p,nstable,npi,nep,nem,ngamma,ptot); - if(p.pdgId()==221) { + if(p.pid()==221) { if(nstable==3 && nem==1 && nem==1 && ngamma==1) { double q = ptot.mass(); double beta = sqrt(1.-4*sqr(me/q)); double p = 1.-sqr(q/meta); double fact = beta*MeV/q*(1.+2.*sqr(me/q))*pow(p,3); - _h_eta->fill(q/MeV,event.weight()/fact); + _h_eta->fill(q/MeV,1./fact); } else if(nstable==2 && ngamma==2) { - _weight_eta += event.weight(); + _weight_eta->fill(); } } else { if(nstable==3 && nem==1 && nem==1 && npi==1) { double q = ptot.mass(); double beta = sqrt(1.-4*sqr(me/q)); double p = sqrt(sqr(1.+sqr(q)/(sqr(momega)-sqr(mpi)))-4.*sqr(momega*q/(sqr(momega)-sqr(mpi)))); double fact = beta*MeV/q*(1.+2.*sqr(me/q))*pow(p,3); - _h_omega->fill(q/MeV,event.weight()/fact); + _h_omega->fill(q/MeV,1./fact); } else if(nstable==2 && ngamma ==1 && npi==1) { - _weight_omega += event.weight(); + _weight_omega->fill(); } } } } /// Normalise histograms etc., after the run void finalize() { static double alpha= 7.2973525664e-3; - scale(_h_eta , 0.75*M_PI/_weight_eta /alpha); - scale(_h_omega, 1.5 *M_PI/_weight_omega/alpha); + scale(_h_eta , 0.75*M_PI/alpha/ *_weight_eta ); + scale(_h_omega, 1.5 *M_PI/alpha/ *_weight_omega); } //@} /// @name Histograms //@{ Histo1DPtr _h_eta,_h_omega; - double _weight_eta,_weight_omega; + CounterPtr _weight_eta,_weight_omega; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(A2_2017_I1486671); } diff --git a/analyses/pluginMisc/E756_2000_I530367.cc b/analyses/pluginMisc/E756_2000_I530367.cc --- a/analyses/pluginMisc/E756_2000_I530367.cc +++ b/analyses/pluginMisc/E756_2000_I530367.cc @@ -1,132 +1,134 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Xi-> Lambda pi asymmetry class E756_2000_I530367 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(E756_2000_I530367); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS" ); // Book histograms - _h_cthetaP = bookHisto1D("cthetaP",20,-1,1); - _h_cthetaM = bookHisto1D("cthetaM",20,-1,1); + book(_h_cthetaP, "cthetaP",20,-1,1); + book(_h_cthetaM, "cthetaM",20,-1,1); } /// Perform the per-event analysis void analyze(const Event& event) { // loop over Xi- baryons - foreach(const Particle& Xi, apply(event, "UFS").particles(Cuts::abspid==3312)) { - int sign = Xi.pdgId()/3312; + for (const Particle& Xi : apply(event, "UFS").particles(Cuts::abspid==3312)) { + int sign = Xi.pid()/3312; if(Xi.children().size()!=2) continue; Particle Lambda,pion1; - if(Xi.children()[0].pdgId()==sign*3122 && - Xi.children()[1].pdgId()==-sign*211) { + if(Xi.children()[0].pid()==sign*3122 && + Xi.children()[1].pid()==-sign*211) { Lambda = Xi.children()[0]; pion1 = Xi.children()[1]; } - else if(Xi.children()[1].pdgId()==sign*3122 && - Xi.children()[0].pdgId()==-sign*211) { + else if(Xi.children()[1].pid()==sign*3122 && + Xi.children()[0].pid()==-sign*211) { Lambda = Xi.children()[1]; pion1 = Xi.children()[0]; } else continue; if(Lambda.children().size()!=2) continue; Particle proton,pion2; - if(Lambda.children()[0].pdgId()==sign*2212 && - Lambda.children()[1].pdgId()==-sign*211) { + if(Lambda.children()[0].pid()==sign*2212 && + Lambda.children()[1].pid()==-sign*211) { proton = Lambda.children()[0]; pion2 = Lambda.children()[1]; } - else if(Lambda.children()[1].pdgId()==sign*2212 && - Lambda.children()[0].pdgId()==-sign*211) { + else if(Lambda.children()[1].pid()==sign*2212 && + Lambda.children()[0].pid()==-sign*211) { proton = Lambda.children()[1]; pion2 = Lambda.children()[0]; } else continue; // boost to xi rest frame first LorentzTransform boost1 = LorentzTransform::mkFrameTransformFromBeta(Xi.momentum().betaVec()); FourMomentum pLambda = boost1.transform(Lambda.momentum()); FourMomentum pproton = boost1.transform(proton.momentum()); // to lambda rest frame LorentzTransform boost2 = LorentzTransform::mkFrameTransformFromBeta(pLambda.betaVec()); Vector3 axis = pLambda.p3().unit(); FourMomentum pp = boost2.transform(pproton); // calculate angle double cTheta = pp.p3().unit().dot(axis); if(sign==1) { - _h_cthetaM->fill(cTheta,1.); + _h_cthetaM->fill(cTheta); } else { - _h_cthetaP->fill(cTheta,1.); + _h_cthetaP->fill(cTheta); } } } pair calcAlpha(Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,0.); double sum1(0.),sum2(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double ai = 0.5*(bin.xMax()-bin.xMin()); double bi = 0.5*ai*(bin.xMax()+bin.xMin()); double Ei = bin.areaErr(); sum1 += sqr(bi/Ei); sum2 += bi/sqr(Ei)*(Oi-ai); } return make_pair(sum2/sum1,sqrt(1./sum1)); } /// Normalise histograms etc., after the run void finalize() { normalize(_h_cthetaP); normalize(_h_cthetaM); // calculate the values of alpha // xibar+ - Scatter2DPtr _h_alphaP = bookScatter2D(1,1,2); + Scatter2DPtr _h_alphaP; + book(_h_alphaP, 1,1,2); pair alpha = calcAlpha(_h_cthetaP); _h_alphaP->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second,alpha.second) ); // xi- - Scatter2DPtr _h_alphaM = bookScatter2D(1,1,1); + Scatter2DPtr _h_alphaM; + book(_h_alphaM, 1,1,1); alpha = calcAlpha(_h_cthetaM); _h_alphaM->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second,alpha.second) ); } //@} /// @name Histograms //@{ Histo1DPtr _h_cthetaP,_h_cthetaM; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(E756_2000_I530367); } diff --git a/analyses/pluginMisc/FOCUS_2006_I693639.cc b/analyses/pluginMisc/FOCUS_2006_I693639.cc --- a/analyses/pluginMisc/FOCUS_2006_I693639.cc +++ b/analyses/pluginMisc/FOCUS_2006_I693639.cc @@ -1,113 +1,114 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Lambda_c -> Lambda pi asymmetry class FOCUS_2006_I693639 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(FOCUS_2006_I693639); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS" ); // Book histograms - _h_ctheta = bookHisto1D("ctheta" , 20,-1,1); + book(_h_ctheta,"ctheta" , 20,-1,1); } /// Perform the per-event analysis void analyze(const Event& event) { // loop over Omega baryons - foreach(const Particle& Lambdac, apply(event, "UFS").particles(Cuts::abspid==4122)) { - int sign = Lambdac.pdgId()/4122; + for( const Particle& Lambdac : apply(event, "UFS").particles(Cuts::abspid==4122)) { + int sign = Lambdac.pid()/4122; if(Lambdac.children().size()!=2) continue; Particle baryon1,meson1; - if(Lambdac.children()[0].pdgId()==sign*3122 && - Lambdac.children()[1].pdgId()==sign*211) { + if(Lambdac.children()[0].pid()==sign*3122 && + Lambdac.children()[1].pid()==sign*211) { baryon1 = Lambdac.children()[0]; meson1 = Lambdac.children()[1]; } - else if(Lambdac.children()[1].pdgId()==sign*3122 && - Lambdac.children()[0].pdgId()==sign*211) { + else if(Lambdac.children()[1].pid()==sign*3122 && + Lambdac.children()[0].pid()==sign*211) { baryon1 = Lambdac.children()[1]; meson1 = Lambdac.children()[0]; } else continue; Particle baryon2,meson2; - if(baryon1.children()[0].pdgId()== sign*2212 && - baryon1.children()[1].pdgId()==-sign*211) { + if(baryon1.children()[0].pid()== sign*2212 && + baryon1.children()[1].pid()==-sign*211) { baryon2 = baryon1.children()[0]; meson2 = baryon1.children()[1]; } - else if(baryon1.children()[1].pdgId()== sign*2212 && - baryon1.children()[0].pdgId()==-sign*211) { + else if(baryon1.children()[1].pid()== sign*2212 && + baryon1.children()[0].pid()==-sign*211) { baryon2 = baryon1.children()[1]; meson2 = baryon1.children()[0]; } else continue; // first boost to the Lambdac rest frame LorentzTransform boost1 = LorentzTransform::mkFrameTransformFromBeta(Lambdac.momentum().betaVec()); FourMomentum pbaryon1 = boost1.transform(baryon1.momentum()); FourMomentum pbaryon2 = boost1.transform(baryon2.momentum()); // to lambda rest frame LorentzTransform boost2 = LorentzTransform::mkFrameTransformFromBeta(pbaryon1.betaVec()); Vector3 axis = pbaryon1.p3().unit(); FourMomentum pp = boost2.transform(pbaryon2); // calculate angle double cTheta = pp.p3().unit().dot(axis); - _h_ctheta->fill(cTheta,1.); + _h_ctheta->fill(cTheta); } } pair calcAlpha(Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,0.); double sum1(0.),sum2(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double ai = 0.5*(bin.xMax()-bin.xMin()); double bi = 0.5*ai*(bin.xMax()+bin.xMin()); double Ei = bin.areaErr(); sum1 += sqr(bi/Ei); sum2 += bi/sqr(Ei)*(Oi-ai); } return make_pair(sum2/sum1,sqrt(1./sum1)); } /// Normalise histograms etc., after the run void finalize() { normalize(_h_ctheta); - Scatter2DPtr _h_alpha = bookScatter2D(1,1,1); + Scatter2DPtr _h_alpha; + book(_h_alpha,1,1,1); pair alpha = calcAlpha(_h_ctheta); _h_alpha->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second,alpha.second) ); } //@} /// @name Histograms //@{ Histo1DPtr _h_ctheta; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(FOCUS_2006_I693639); } diff --git a/analyses/pluginMisc/HyperCP_2005_I677384.cc b/analyses/pluginMisc/HyperCP_2005_I677384.cc --- a/analyses/pluginMisc/HyperCP_2005_I677384.cc +++ b/analyses/pluginMisc/HyperCP_2005_I677384.cc @@ -1,134 +1,137 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief asymmetry in Omega-> Lambda K class HyperCP_2005_I677384 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(HyperCP_2005_I677384); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS" ); // Book histograms - _h_cthetaP = bookHisto1D("cthetaP" ,20,-1,1); - _h_cthetaM = bookHisto1D("cthetaM" ,20,-1,1); - _h_cthetaAll = bookHisto1D("cthetaAll",20,-1,1); + book(_h_cthetaP , "cthetaP" ,20,-1,1); + book(_h_cthetaM , "cthetaM" ,20,-1,1); + book(_h_cthetaAll, "cthetaAll",20,-1,1); } /// Perform the per-event analysis void analyze(const Event& event) { // loop over Omega baryons - foreach(const Particle& Omega, apply(event, "UFS").particles(Cuts::abspid==3334)) { - int sign = Omega.pdgId()/3334; + for(const Particle& Omega : apply(event, "UFS").particles(Cuts::abspid==3334)) { + int sign = Omega.pid()/3334; if(Omega.children().size()!=2) continue; Particle Lambda,kaon; - if(Omega.children()[0].pdgId()==sign*3122 && - Omega.children()[1].pdgId()==-sign*321) { + if(Omega.children()[0].pid()==sign*3122 && + Omega.children()[1].pid()==-sign*321) { Lambda = Omega.children()[0]; kaon = Omega.children()[1]; } - else if(Omega.children()[1].pdgId()==sign*3122 && - Omega.children()[0].pdgId()==-sign*321) { + else if(Omega.children()[1].pid()==sign*3122 && + Omega.children()[0].pid()==-sign*321) { Lambda = Omega.children()[1]; kaon = Omega.children()[0]; } else continue; if(Lambda.children().size()!=2) continue; Particle proton,pion; - if(Lambda.children()[0].pdgId()==sign*2212 && - Lambda.children()[1].pdgId()==-sign*211) { + if(Lambda.children()[0].pid()==sign*2212 && + Lambda.children()[1].pid()==-sign*211) { proton = Lambda.children()[0]; pion = Lambda.children()[1]; } - else if(Lambda.children()[1].pdgId()==sign*2212 && - Lambda.children()[0].pdgId()==-sign*211) { + else if(Lambda.children()[1].pid()==sign*2212 && + Lambda.children()[0].pid()==-sign*211) { proton = Lambda.children()[1]; pion = Lambda.children()[0]; } else continue; // first boost to the Omega rest frame LorentzTransform boost1 = LorentzTransform::mkFrameTransformFromBeta(Omega.momentum().betaVec()); FourMomentum pLambda = boost1.transform(Lambda.momentum()); FourMomentum pproton = boost1.transform(proton.momentum()); // to lambda rest frame LorentzTransform boost2 = LorentzTransform::mkFrameTransformFromBeta(pLambda.betaVec()); Vector3 axis = pLambda.p3().unit(); FourMomentum pp = boost2.transform(pproton); // calculate angle double cTheta = pp.p3().unit().dot(axis); _h_cthetaAll->fill(cTheta,1.); if(sign==1) { _h_cthetaM->fill(cTheta,1.); } else { _h_cthetaP->fill(cTheta,1.); } } } pair calcAlpha(Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,0.); double sum1(0.),sum2(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double ai = 0.5*(bin.xMax()-bin.xMin()); double bi = 0.5*ai*(bin.xMax()+bin.xMin()); double Ei = bin.areaErr(); sum1 += sqr(bi/Ei); sum2 += bi/sqr(Ei)*(Oi-ai); } return make_pair(sum2/sum1,sqrt(1./sum1)); } /// Normalise histograms etc., after the run void finalize() { normalize(_h_cthetaP ); normalize(_h_cthetaM ); normalize(_h_cthetaAll); // calculate the values of alpha - Scatter2DPtr _h_alphaP = bookScatter2D(1,1,1); + Scatter2DPtr _h_alphaP; + book(_h_alphaP,1,1,1); pair alpha = calcAlpha(_h_cthetaP); _h_alphaP->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second,alpha.second) ); - Scatter2DPtr _h_alphaM = bookScatter2D(1,1,2); + Scatter2DPtr _h_alphaM; + book(_h_alphaM,1,1,2); alpha = calcAlpha(_h_cthetaM); _h_alphaM->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second,alpha.second) ); - Scatter2DPtr _h_alphaAll = bookScatter2D(1,1,3); + Scatter2DPtr _h_alphaAll; + book(_h_alphaAll,1,1,3); alpha = calcAlpha(_h_cthetaAll); _h_alphaAll->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second,alpha.second) ); } //@} /// @name Histograms //@{ Histo1DPtr _h_cthetaP,_h_cthetaM,_h_cthetaAll; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(HyperCP_2005_I677384); } diff --git a/analyses/pluginMisc/NA48_2010_I868871.cc b/analyses/pluginMisc/NA48_2010_I868871.cc --- a/analyses/pluginMisc/NA48_2010_I868871.cc +++ b/analyses/pluginMisc/NA48_2010_I868871.cc @@ -1,226 +1,231 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" - +#include namespace Rivet { /// @brief asymmetrics in Xi0 decays class NA48_2010_I868871 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(NA48_2010_I868871); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS" ); // Book histograms - _h_ctheta_pi0 = bookHisto1D("ctheta_pi0" , 20,-1,1); - _h_ctheta_gamma = bookHisto1D("ctheta_gamma", 20,-1,1); + book(_h_ctheta_pi0 , "ctheta_pi0" , 20,-1,1); + book(_h_ctheta_gamma, "ctheta_gamma", 20,-1,1); double step=0.1; double xmin=-1.; for(unsigned int ix=0;ix<20;++ix) { - ostringstream title; + Histo1DPtr temp; + std::ostringstream title; title << "ctheta_Sigma_" << ix; - _h_ctheta_Sigma.addHistogram(xmin, xmin+step, bookHisto1D(title.str(), 20,-1,1)); + book(temp,title.str(), 20,-1,1); + _h_ctheta_Sigma.add(xmin, xmin+step, temp); xmin+=step; } _nSigma=0.; } /// Perform the per-event analysis void analyze(const Event& event) { // loop over Omega baryons - foreach(const Particle& Xi, apply(event, "UFS").particles(Cuts::abspid==3322)) { - int sign = Xi.pdgId()/3322; + for(const Particle& Xi : apply(event, "UFS").particles(Cuts::abspid==3322)) { + int sign = Xi.pid()/3322; if(Xi.children().size()!=2) continue; Particle baryon1,meson1; unsigned int mode(0); - if(Xi.children()[0].pdgId()==sign*3122 && - Xi.children()[1].pdgId()==111) { + if(Xi.children()[0].pid()==sign*3122 && + Xi.children()[1].pid()==111) { baryon1 = Xi.children()[0]; meson1 = Xi.children()[1]; mode=1; } - else if(Xi.children()[1].pdgId()==sign*3122 && - Xi.children()[0].pdgId()==111) { + else if(Xi.children()[1].pid()==sign*3122 && + Xi.children()[0].pid()==111) { baryon1 = Xi.children()[1]; meson1 = Xi.children()[0]; mode=1; } - else if(Xi.children()[0].pdgId()==sign*3122 && - Xi.children()[1].pdgId()==22) { + else if(Xi.children()[0].pid()==sign*3122 && + Xi.children()[1].pid()==22) { baryon1 = Xi.children()[0]; meson1 = Xi.children()[1]; mode=2; } - else if(Xi.children()[1].pdgId()==sign*3122 && - Xi.children()[0].pdgId()==22) { + else if(Xi.children()[1].pid()==sign*3122 && + Xi.children()[0].pid()==22) { baryon1 = Xi.children()[1]; meson1 = Xi.children()[0]; mode=2; } - else if(Xi.children()[0].pdgId()==sign*3212 && - Xi.children()[1].pdgId()==22) { + else if(Xi.children()[0].pid()==sign*3212 && + Xi.children()[1].pid()==22) { baryon1 = Xi.children()[0]; meson1 = Xi.children()[1]; mode=3; } - else if(Xi.children()[1].pdgId()==sign*3212 && - Xi.children()[0].pdgId()==22) { + else if(Xi.children()[1].pid()==sign*3212 && + Xi.children()[0].pid()==22) { baryon1 = Xi.children()[1]; meson1 = Xi.children()[0]; mode=3; } else continue; if(baryon1.children().size()!=2) continue; Particle baryon2,meson2,baryon3,meson3; if(mode==1 || mode ==2) { - if(baryon1.children()[0].pdgId()== sign*2212 && - baryon1.children()[1].pdgId()==-sign*211) { + if(baryon1.children()[0].pid()== sign*2212 && + baryon1.children()[1].pid()==-sign*211) { baryon2 = baryon1.children()[0]; meson2 = baryon1.children()[1]; } - else if(baryon1.children()[1].pdgId()== sign*2212 && - baryon1.children()[0].pdgId()==-sign*211) { + else if(baryon1.children()[1].pid()== sign*2212 && + baryon1.children()[0].pid()==-sign*211) { baryon2 = baryon1.children()[1]; meson2 = baryon1.children()[0]; } else continue; } else if(mode==3) { - if(baryon1.children()[0].pdgId()== sign*3122 && - baryon1.children()[1].pdgId()== 22) { + if(baryon1.children()[0].pid()== sign*3122 && + baryon1.children()[1].pid()== 22) { baryon2 = baryon1.children()[0]; meson2 = baryon1.children()[1]; } - else if(baryon1.children()[1].pdgId()== sign*3122 && - baryon1.children()[0].pdgId()== 22) { + else if(baryon1.children()[1].pid()== sign*3122 && + baryon1.children()[0].pid()== 22) { baryon2 = baryon1.children()[1]; meson2 = baryon1.children()[0]; } else continue; - if(baryon2.children()[0].pdgId()== sign*2212 && - baryon2.children()[1].pdgId()==-sign*211) { + if(baryon2.children()[0].pid()== sign*2212 && + baryon2.children()[1].pid()==-sign*211) { baryon3 = baryon2.children()[0]; meson3 = baryon2.children()[1]; } - else if(baryon2.children()[1].pdgId()== sign*2212 && - baryon2.children()[0].pdgId()==-sign*211) { + else if(baryon2.children()[1].pid()== sign*2212 && + baryon2.children()[0].pid()==-sign*211) { baryon3 = baryon2.children()[1]; meson3 = baryon2.children()[0]; } else continue; } // first boost to the Xi rest frame LorentzTransform boost1 = LorentzTransform::mkFrameTransformFromBeta(Xi.momentum().betaVec()); FourMomentum pbaryon1 = boost1.transform(baryon1.momentum()); FourMomentum pbaryon2 = boost1.transform(baryon2.momentum()); // to lambda rest frame LorentzTransform boost2 = LorentzTransform::mkFrameTransformFromBeta(pbaryon1.betaVec()); Vector3 axis = pbaryon1.p3().unit(); FourMomentum pp = boost2.transform(pbaryon2); // calculate angle double cTheta = pp.p3().unit().dot(axis); if(mode==1) { _h_ctheta_pi0->fill(cTheta,1.); } else if(mode==2) { _h_ctheta_gamma->fill(cTheta,1.); } else if(mode==3) { FourMomentum pbaryon3 = boost1.transform(baryon3.momentum()); FourMomentum pp2 = boost2.transform(pbaryon3); Vector3 axis2 = pp.p3().unit(); double cTheta2 = pp2.p3().unit().dot(axis2); _h_ctheta_Sigma.fill(cTheta,cTheta2,1.); _nSigma += 1.; } } } pair calcAlpha(Histo1DPtr hist) { if(hist->numEntries()==0.) return make_pair(0.,0.); double sum1(0.),sum2(0.); for (auto bin : hist->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double ai = 0.5*(bin.xMax()-bin.xMin()); double bi = 0.5*ai*(bin.xMax()+bin.xMin()); double Ei = bin.areaErr(); sum1 += sqr(bi/Ei); sum2 += bi/sqr(Ei)*(Oi-ai); } return make_pair(sum2/sum1,sqrt(1./sum1)); } - pair calcAlpha(BinnedHistogram hist) { + pair calcAlpha(BinnedHistogram hist) { double sum1(0.),sum2(0.); double step=0.1; double xmin=-1.; for(unsigned int ix=0;ix<20;++ix) { double xsum=2.*xmin+step; Histo1DPtr h2 = hist.histo(xmin+0.5*step); for (auto bin : h2->bins() ) { double Oi = bin.area(); if(Oi==0.) continue; double ai = 0.25*(bin.xMax()-bin.xMin())*step; double bi = 0.25*ai*(bin.xMax()+bin.xMin())*xsum; double Ei = bin.areaErr(); sum1 += sqr(bi/Ei); sum2 += bi/sqr(Ei)*(Oi-ai); } xmin+=step; } return make_pair(sum2/sum1,sqrt(1./sum1)); } /// Normalise histograms etc., after the run void finalize() { // Xi0 -> Lambda0 pi0 normalize(_h_ctheta_pi0); - Scatter2DPtr _h_alpha_pi0 = bookScatter2D(1,1,1); + Scatter2DPtr _h_alpha_pi0; + book(_h_alpha_pi0,1,1,1); pair alpha = calcAlpha(_h_ctheta_pi0); _h_alpha_pi0->addPoint(0.5, alpha.first, make_pair(0.5,0.5), make_pair(alpha.second,alpha.second) ); // Xi0 -> Lambda gamma (N.B. sign due defns) normalize(_h_ctheta_gamma); - Scatter2DPtr _h_alpha_gamma = bookScatter2D(2,1,1); + Scatter2DPtr _h_alpha_gamma; + book(_h_alpha_gamma,2,1,1); alpha = calcAlpha(_h_ctheta_gamma); _h_alpha_gamma->addPoint(0.5,-alpha.first, make_pair(0.5,0.5), make_pair(alpha.second,alpha.second) ); // Xi0 -> Sigma gamma _h_ctheta_Sigma.scale(1./_nSigma,this); - Scatter2DPtr _h_alpha_Sigma = bookScatter2D(3,1,1); + Scatter2DPtr _h_alpha_Sigma; + book(_h_alpha_Sigma,3,1,1); alpha = calcAlpha(_h_ctheta_Sigma); _h_alpha_Sigma->addPoint(0.5,alpha.first, make_pair(0.5,0.5), make_pair(alpha.second,alpha.second) ); } //@} /// @name Histograms //@{ Histo1DPtr _h_ctheta_pi0,_h_ctheta_gamma; - BinnedHistogram _h_ctheta_Sigma; + BinnedHistogram _h_ctheta_Sigma; double _nSigma; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(NA48_2010_I868871); } diff --git a/analyses/pluginNovosibirsk/CMD2_2000_I532970.cc b/analyses/pluginNovosibirsk/CMD2_2000_I532970.cc --- a/analyses/pluginNovosibirsk/CMD2_2000_I532970.cc +++ b/analyses/pluginNovosibirsk/CMD2_2000_I532970.cc @@ -1,166 +1,166 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class CMD2_2000_I532970 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(CMD2_2000_I532970); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); book(_numOmegaPiPi, "TMP/OmegaPiPi"); book(_numEtaPiPi, "TMP/EtaPiPi"); book(_num5Pi, "TMP/5Pi"); } void findChildren(const Particle & p,map & nRes, int &ncount) { for (const Particle &child : p.children()) { if(child.children().empty()) { --nRes[child.pid()]; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); for (const Particle& p : fs.particles()) { nCount[p.pid()] += 1; ++ntotal; } const FinalState& ufs = apply(event, "UFS"); bool foundRes = false; for (const Particle& p : ufs.particles()) { if(p.children().empty()) continue; // find the eta if(p.pid()==221) { map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); // eta pi+pi- if(ncount!=2) continue; bool matched = true; for(auto const & val : nRes) { if(abs(val.first)==211) { if(val.second !=1) { matched = false; break; } } else if(val.second!=0) { matched = false; break; } } if(matched) { _numEtaPiPi->fill(); foundRes = true; } } // find the omega - else if(p.pdgId()==223) { + else if(p.pid()==223) { map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); // eta pi+pi- if(ncount!=2) continue; bool matched = true; for(auto const & val : nRes) { if(abs(val.first)==211) { if(val.second !=1) { matched = false; break; } } else if(val.second!=0) { matched = false; break; } } if(matched) { _numOmegaPiPi->fill(); foundRes = true; } } } if(foundRes) vetoEvent; if(nCount[-211]==2&&nCount[211]==2&&nCount[111]==1) _num5Pi->fill(); } /// Normalise histograms etc., after the run void finalize() { for(unsigned int ix=1;ix<4;++ix) { double sigma,error; if(ix==1) { sigma = _numOmegaPiPi->val(); error = _numOmegaPiPi->err(); } else if(ix==2) { sigma = _numEtaPiPi->val(); error = _numEtaPiPi->err(); } else { sigma = _num5Pi->val(); error = _num5Pi->err(); } sigma *= crossSection()/ sumOfWeights() /nanobarn; error *= crossSection()/ sumOfWeights() /nanobarn; Scatter2D temphisto(refData(ix, 1, 1)); Scatter2DPtr mult; book(mult, ix, 1, 1); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/GeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } } //@} /// @name Histograms //@{ CounterPtr _numEtaPiPi,_numOmegaPiPi,_num5Pi; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(CMD2_2000_I532970); } diff --git a/analyses/pluginNovosibirsk/CMD3_2013_I1217420.cc b/analyses/pluginNovosibirsk/CMD3_2013_I1217420.cc --- a/analyses/pluginNovosibirsk/CMD3_2013_I1217420.cc +++ b/analyses/pluginNovosibirsk/CMD3_2013_I1217420.cc @@ -1,100 +1,100 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class CMD3_2013_I1217420 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(CMD3_2013_I1217420); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); // Book histograms - _c_all = bookCounter("/TMP/all"); + book(_c_all, "/TMP/all"); } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for(const Particle &child : p.children()) { if(child.children().empty()) { - --nRes[child.pdgId()]; + --nRes[child.pid()]; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { // find the final-state particles - double weight = event.weight(); const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } if(ntotal==6 && nCount[211]==3 && nCount[-211]==3) { - _c_all->fill(weight); + _c_all->fill(); } } /// Normalise histograms etc., after the run void finalize() { double fact = crossSection()/nanobarn/sumOfWeights(); double sigma = _c_all->val()*fact; double error = _c_all->err()*fact; Scatter2D temphisto(refData(1, 1, 1)); - Scatter2DPtr mult = bookScatter2D(1, 1, 1); + Scatter2DPtr mult; + book(mult, 1, 1, 1); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/MeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } //@} /// @name Histograms //@{ CounterPtr _c_all; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(CMD3_2013_I1217420); } diff --git a/analyses/pluginNovosibirsk/CMD3_2017_I1606078.cc b/analyses/pluginNovosibirsk/CMD3_2017_I1606078.cc --- a/analyses/pluginNovosibirsk/CMD3_2017_I1606078.cc +++ b/analyses/pluginNovosibirsk/CMD3_2017_I1606078.cc @@ -1,204 +1,204 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief CMD3 pi+pi-pi0eta cross section class CMD3_2017_I1606078 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(CMD3_2017_I1606078); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); - _c_all = bookCounter("/TMP/all"); - _c_omega = bookCounter("/TMP/omega"); - _c_rho = bookCounter("/TMP/rho"); - _c_other = bookCounter("/TMP/other"); + book(_c_all , "/TMP/all"); + book(_c_omega, "/TMP/omega"); + book(_c_rho , "/TMP/rho"); + book(_c_other, "/TMP/other"); } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for(const Particle &child : p.children()) { if(child.children().empty()) { - --nRes[child.pdgId()]; + --nRes[child.pid()]; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { // find the final-state particles - double weight = event.weight(); const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } // find omega/phi + eta const FinalState& ufs = apply(event, "UFS"); bool found = false, foundOmegaPhi = false; - foreach (const Particle& p, ufs.particles()) { + for (const Particle& p : ufs.particles()) { if(p.children().empty()) continue; // find the eta - if(p.pdgId()!=221) continue; + if(p.pid()!=221) continue; map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); // eta pi+pi-pi0 if(ncount==3) { bool matched = true; for(auto const & val : nRes) { if(abs(val.first)==211 || val.first==111 ) { if(val.second !=1) { matched = false; break; } } else if(val.second!=0) { matched = false; break; } } if(matched) { - _c_all->fill(event.weight()); + _c_all->fill(); found = true; } } - foreach (const Particle& p2, ufs.particles()) { - if(p2.pdgId()!=223 && p2.pdgId()!=333) continue; + for (const Particle& p2 : ufs.particles()) { + if(p2.pid()!=223 && p2.pid()!=333) continue; map nResB = nRes; int ncountB = ncount; findChildren(p2,nResB,ncountB); if(ncountB!=0) continue; bool matched2 = true; for(auto const & val : nResB) { if(val.second!=0) { matched2 = false; break; } } if(matched2) { - if(p2.pdgId()==223) - _c_omega->fill(event.weight()); + if(p2.pid()==223) + _c_omega->fill(); foundOmegaPhi=true; } } } // find a_0 rho bool founda0Rho=false; - foreach (const Particle& p, ufs.particles()) { + for (const Particle& p : ufs.particles()) { if(p.children().empty()) continue; // find the rho - if(p.pdgId()!=113 && abs(p.pdgId())!=213) + if(p.pid()!=113 && abs(p.pid())!=213) continue; map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); int a1id(0); - if(p.pdgId()==213) + if(p.pid()==213) a1id = 9000211; - else if(p.pdgId()==-213) + else if(p.pid()==-213) a1id =-9000211; else a1id = 9000111; - foreach (const Particle& p2, ufs.particles()) { - if(p2.pdgId()!=a1id) continue; + for (const Particle& p2 : ufs.particles()) { + if(p2.pid()!=a1id) continue; map nResB = nRes; int ncountB = ncount; findChildren(p2,nResB,ncountB); if(ncountB!=0) continue; bool matched = true; for(auto const & val : nResB) { if(val.second!=0) { matched = false; break; } } if(matched) { - _c_rho->fill(event.weight()); + _c_rho->fill(); founda0Rho=true; break; } } if(founda0Rho) break; } if(found && !foundOmegaPhi && ! founda0Rho) - _c_other->fill(weight); + _c_other->fill(); } /// Normalise histograms etc., after the run void finalize() { double fact = crossSection()/nanobarn/sumOfWeights(); for(unsigned int ix=1;ix<5;++ix) { double sigma(0.),error(0.); if(ix==1) { sigma = _c_all->val()*fact; error = _c_all->err()*fact; } else if(ix==2) { sigma = _c_omega->val()*fact; error = _c_omega->err()*fact; } else if(ix==3) { sigma = _c_rho->val()*fact; error = _c_rho->err()*fact; } else if(ix==4) { sigma = _c_other->val()*fact; error = _c_other->err()*fact; } Scatter2D temphisto(refData(1, 1, ix)); - Scatter2DPtr mult = bookScatter2D(1, 1, ix); + Scatter2DPtr mult; + book(mult, 1, 1, ix); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/MeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } } //@} /// @name Histograms //@{ CounterPtr _c_all, _c_omega, _c_rho, _c_other; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(CMD3_2017_I1606078); } diff --git a/analyses/pluginNovosibirsk/CMD3_2019_I1720610.cc b/analyses/pluginNovosibirsk/CMD3_2019_I1720610.cc --- a/analyses/pluginNovosibirsk/CMD3_2019_I1720610.cc +++ b/analyses/pluginNovosibirsk/CMD3_2019_I1720610.cc @@ -1,151 +1,151 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class CMD3_2019_I1720610 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(CMD3_2019_I1720610); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); // Book histograms - _c_all = bookCounter("/TMP/all"); - _c_omega = bookCounter("/TMP/omega"); - _c_eta = bookCounter("/TMP/eta"); + book(_c_all , "/TMP/all"); + book(_c_omega, "/TMP/omega"); + book(_c_eta , "/TMP/eta"); } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for(const Particle &child : p.children()) { if(child.children().empty()) { - --nRes[child.pdgId()]; + --nRes[child.pid()]; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { // find the final-state particles - double weight = event.weight(); const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } if(ntotal==7 && nCount[211]==3 && nCount[-211]==3 && nCount[111] ==1 ) { - _c_all->fill(weight); + _c_all->fill(); } // find omega/phi + eta const FinalState& ufs = apply(event, "UFS"); bool found=false; - foreach (const Particle& p, ufs.particles()) { + for (const Particle& p : ufs.particles()) { if(p.children().empty()) continue; // find the eta/omega - if(p.pdgId()!=221 && p.pdgId()!=223) continue; + if(p.pid()!=221 && p.pid()!=223) continue; map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); // eta/omega 2(pi+pi-) if(ncount==4) { bool matched = true; for(auto const & val : nRes) { if(abs(val.first)==211 ) { if(val.second !=2) { matched = false; break; } } else if(val.second!=0) { matched = false; break; } } if(matched) { - if(p.pdgId()==221) - _c_eta->fill(event.weight()); - else if(p.pdgId()==223) - _c_omega->fill(event.weight()); + if(p.pid()==221) + _c_eta->fill(); + else if(p.pid()==223) + _c_omega->fill(); found = true; break; } } if(found) break; } } /// Normalise histograms etc., after the run void finalize() { double fact = crossSection()/nanobarn/sumOfWeights(); for(unsigned int ix=1;ix<4;++ix) { double sigma(0.),error(0.); if(ix==1) { sigma = _c_all->val()*fact; error = _c_all->err()*fact; } else if(ix==2) { sigma = _c_eta->val()*fact; error = _c_eta->err()*fact; } else if(ix==3) { sigma = _c_omega->val()*fact; error = _c_omega->err()*fact; } Scatter2D temphisto(refData(1, 1, ix)); - Scatter2DPtr mult = bookScatter2D(1, 1, ix); + Scatter2DPtr mult; + book(mult, 1, 1, ix); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/MeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } } //@} /// @name Histograms //@{ CounterPtr _c_all,_c_omega,_c_eta; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(CMD3_2019_I1720610); } diff --git a/analyses/pluginNovosibirsk/SND_2000_I525398.cc b/analyses/pluginNovosibirsk/SND_2000_I525398.cc --- a/analyses/pluginNovosibirsk/SND_2000_I525398.cc +++ b/analyses/pluginNovosibirsk/SND_2000_I525398.cc @@ -1,97 +1,96 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief phi -> eta pi0 gamma class SND_2000_I525398 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(SND_2000_I525398); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { declare(UnstableParticles(), "UFS"); - _h_etapi = bookHisto1D( 1, 1, 1); - _nPhi=0.; + book(_h_etapi, 1, 1, 1); + book(_nPhi,"TMP/nPhi"); } void findDecayProducts(const Particle & mother, unsigned int & nstable, unsigned int & neta, unsigned int & npi, unsigned int & ngamma, FourMomentum & ptot) { for(const Particle & p : mother.children()) { - int id = p.pdgId(); + int id = p.pid(); if ( id == PID::ETA ) { ++neta; ++nstable; ptot += p.momentum(); } else if (id == PID::PI0) { ++npi; ++nstable; ptot += p.momentum(); } else if (id == PID::GAMMA) { ++ngamma; ++nstable; } else if (id == PID::PIPLUS || id == PID::PIMINUS) { ++nstable; } else if ( !p.children().empty() ) { findDecayProducts(p, nstable, neta, npi, ngamma, ptot); } else ++nstable; } } /// Perform the per-event analysis void analyze(const Event& event) { // Loop over phis for(const Particle& phi : apply(event, "UFS").particles(Cuts::abspid==PID::PHI)) { - _nPhi+=event.weight(); + _nPhi->fill(); unsigned int nstable(0),neta(0),npi(0),ngamma(0); FourMomentum p_tot(0,0,0,0); findDecayProducts(phi, nstable, neta, npi, ngamma, p_tot); if(nstable!=3) continue; if(neta==1 && npi==1 && ngamma==1 ) { - _h_etapi->fill(p_tot.mass()/MeV, event.weight()); + _h_etapi->fill(p_tot.mass()/MeV); } } } /// Normalise histograms etc., after the run void finalize() { // normalise to total no of phi mesons - // and mult by 10^7 due normalisation in paper - scale( _h_etapi, 1./_nPhi); + scale( _h_etapi, 1./ *_nPhi); } //@} /// @name Histograms //@{ Histo1DPtr _h_etapi; - double _nPhi; + CounterPtr _nPhi; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(SND_2000_I525398); } diff --git a/analyses/pluginNovosibirsk/SND_2000_I527094.cc b/analyses/pluginNovosibirsk/SND_2000_I527094.cc --- a/analyses/pluginNovosibirsk/SND_2000_I527094.cc +++ b/analyses/pluginNovosibirsk/SND_2000_I527094.cc @@ -1,89 +1,89 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief phi -> pi0 pi0 gamma class SND_2000_I527094 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(SND_2000_I527094); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { declare(UnstableParticles(), "UFS"); - _h_pipi = bookHisto1D( 1, 1, 1); - _nPhi=0.; + book(_h_pipi, 1, 1, 1); + book(_nPhi,"TMP/nPhi"); } void findDecayProducts(const Particle & mother, unsigned int & nstable, unsigned int & npi, unsigned int & ngamma, FourMomentum & ptot) { for(const Particle & p : mother.children()) { - int id = p.pdgId(); + int id = p.pid(); if (id == PID::PI0) { ++npi; ++nstable; ptot += p.momentum(); } else if (id == PID::GAMMA) { ++ngamma; ++nstable; } else if (id == PID::PIPLUS || id == PID::PIMINUS) { ++nstable; } else if ( !p.children().empty() ) { findDecayProducts(p, nstable, npi, ngamma, ptot); } else ++nstable; } } /// Perform the per-event analysis void analyze(const Event& event) { // Loop over phis for(const Particle& phi : apply(event, "UFS").particles(Cuts::abspid==PID::PHI)) { - _nPhi+=event.weight(); + _nPhi->fill(); unsigned int nstable(0),npi(0),ngamma(0); FourMomentum p_tot(0,0,0,0); findDecayProducts(phi, nstable, npi, ngamma, p_tot); if(nstable!=3) continue; if(npi==2 && ngamma==1 ) { - _h_pipi->fill(p_tot.mass()/MeV, event.weight()); + _h_pipi->fill(p_tot.mass()/MeV); } } } /// Normalise histograms etc., after the run void finalize() { - scale( _h_pipi, 1./_nPhi); + scale( _h_pipi, 1./ *_nPhi); } //@} /// @name Histograms //@{ Histo1DPtr _h_pipi; - double _nPhi; + CounterPtr _nPhi; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(SND_2000_I527094); } diff --git a/analyses/pluginNovosibirsk/SND_2001_I558279.cc b/analyses/pluginNovosibirsk/SND_2001_I558279.cc --- a/analyses/pluginNovosibirsk/SND_2001_I558279.cc +++ b/analyses/pluginNovosibirsk/SND_2001_I558279.cc @@ -1,93 +1,92 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief phi -> 3 pion decay class SND_2001_I558279 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(SND_2001_I558279); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(UnstableParticles(), "UFS"); // Book histograms - _h_pm = bookHisto1D(1, 1, 1); - _h_p0 = bookHisto1D(2, 1, 1); + book(_h_pm, 1, 1, 1); + book(_h_p0, 2, 1, 1); } void findDecayProducts(const Particle & mother, unsigned int & nstable, Particles &pip, Particles &pim, Particles &pi0) { for(const Particle & p : mother.children()) { - int id = p.pdgId(); + int id = p.pid(); if (id == PID::PIMINUS ) { pim.push_back(p); ++nstable; } else if (id == PID::PIPLUS) { pip.push_back(p); ++nstable; } else if (id == PID::PI0) { pi0.push_back(p); ++nstable; } else if ( !p.children().empty() ) { findDecayProducts(p, nstable, pip,pim,pi0); } else ++nstable; } } /// Perform the per-event analysis void analyze(const Event& event) { - const double weight=event.weight(); - foreach(const Particle& p, apply(event, "UFS").particles(Cuts::pid==333)) { + for (const Particle& p : apply(event, "UFS").particles(Cuts::pid==333)) { Particles pip,pim,pi0; unsigned int nstable(0); findDecayProducts(p, nstable, pip,pim,pi0); if(nstable==3 && pip.size()==1 && pim.size()==1&&pi0.size()==1) { - _h_pm->fill((pip[0].momentum()+pim[0].momentum()).mass()/MeV,weight); - _h_p0->fill((pip[0].momentum()+pi0[0].momentum()).mass()/MeV,weight); - _h_p0->fill((pim[0].momentum()+pi0[0].momentum()).mass()/MeV,weight); + _h_pm->fill((pip[0].momentum()+pim[0].momentum()).mass()/MeV); + _h_p0->fill((pip[0].momentum()+pi0[0].momentum()).mass()/MeV); + _h_p0->fill((pim[0].momentum()+pi0[0].momentum()).mass()/MeV); } } } /// Normalise histograms etc., after the run void finalize() { normalize(_h_pm); normalize(_h_p0); } //@} /// @name Histograms //@{ Histo1DPtr _h_pm,_h_p0; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(SND_2001_I558279); } diff --git a/analyses/pluginNovosibirsk/SND_2015_I1332929.cc b/analyses/pluginNovosibirsk/SND_2015_I1332929.cc --- a/analyses/pluginNovosibirsk/SND_2015_I1332929.cc +++ b/analyses/pluginNovosibirsk/SND_2015_I1332929.cc @@ -1,121 +1,121 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief Add a short analysis description here class SND_2015_I1332929 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(SND_2015_I1332929); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); book(_numEtaPiPi, "TMP/EtaPiPi"); } void findChildren(const Particle & p,map & nRes, int &ncount) { for (const Particle &child : p.children()) { if(child.children().empty()) { --nRes[child.pid()]; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); for (const Particle& p : fs.particles()) { nCount[p.pid()] += 1; ++ntotal; } const FinalState& ufs = apply(event, "UFS"); for (const Particle& p : ufs.particles()) { if(p.children().empty()) continue; // find the eta - if(p.pdgId()==221) { + if(p.pid()==221) { map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); // eta pi+pi- if(ncount!=2) continue; bool matched = true; for(auto const & val : nRes) { if(abs(val.first)==211) { if(val.second !=1) { matched = false; break; } } else if(val.second!=0) { matched = false; break; } } if(matched) _numEtaPiPi->fill(); } } } /// Normalise histograms etc., after the run void finalize() { double sigma = _numEtaPiPi->val(); double error = _numEtaPiPi->err(); sigma *= crossSection()/ sumOfWeights() /nanobarn; error *= crossSection()/ sumOfWeights() /nanobarn; Scatter2D temphisto(refData(1, 1, 1)); Scatter2DPtr mult; book(mult, 1, 1, 1); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/GeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } //@} /// @name Histograms //@{ CounterPtr _numEtaPiPi; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(SND_2015_I1332929); } diff --git a/analyses/pluginNovosibirsk/SND_2019_I1726419.cc b/analyses/pluginNovosibirsk/SND_2019_I1726419.cc --- a/analyses/pluginNovosibirsk/SND_2019_I1726419.cc +++ b/analyses/pluginNovosibirsk/SND_2019_I1726419.cc @@ -1,169 +1,169 @@ // -*- C++ -*- #include "Rivet/Analysis.hh" #include "Rivet/Projections/FinalState.hh" #include "Rivet/Projections/UnstableParticles.hh" namespace Rivet { /// @brief SND pi+pi-pi0eta cross section class SND_2019_I1726419 : public Analysis { public: /// Constructor DEFAULT_RIVET_ANALYSIS_CTOR(SND_2019_I1726419); /// @name Analysis methods //@{ /// Book histograms and initialise projections before the run void init() { // Initialise and register projections declare(FinalState(), "FS"); declare(UnstableParticles(), "UFS"); - _c_all = bookCounter("/TMP/all"); - _c_omega = bookCounter("/TMP/omega"); - _c_phi = bookCounter("/TMP/phi"); - _c_rho = bookCounter("/TMP/rho"); + book(_c_all , "/TMP/all"); + book(_c_omega, "/TMP/omega"); + book(_c_phi , "/TMP/phi"); + book(_c_rho , "/TMP/rho"); } void findChildren(const Particle & p,map & nRes, int &ncount) { - foreach(const Particle &child, p.children()) { + for( const Particle &child : p.children()) { if(child.children().empty()) { - --nRes[child.pdgId()]; + --nRes[child.pid()]; --ncount; } else findChildren(child,nRes,ncount); } } /// Perform the per-event analysis void analyze(const Event& event) { - double weight = event.weight(); const FinalState& fs = apply(event, "FS"); map nCount; int ntotal(0); - foreach (const Particle& p, fs.particles()) { - nCount[p.pdgId()] += 1; + for (const Particle& p : fs.particles()) { + nCount[p.pid()] += 1; ++ntotal; } const FinalState& ufs = apply(event, "UFS"); bool found = false, foundOmegaPhi = false; - foreach (const Particle& p, ufs.particles()) { + for (const Particle& p : ufs.particles()) { if(p.children().empty()) continue; // find the eta - if(p.pdgId()!=221) continue; + if(p.pid()!=221) continue; map nRes = nCount; int ncount = ntotal; findChildren(p,nRes,ncount); // eta pi+pi- if(ncount==3) { bool matched = true; for(auto const & val : nRes) { if(abs(val.first)==211 || val.first==111 ) { if(val.second !=1) { matched = false; break; } } else if(val.second!=0) { matched = false; break; } } if(matched) { - _c_all->fill(event.weight()); + _c_all->fill(); found = true; } } - foreach (const Particle& p2, ufs.particles()) { - if(p2.pdgId()!=223 && p2.pdgId()!=333) continue; + for (const Particle& p2 : ufs.particles()) { + if(p2.pid()!=223 && p2.pid()!=333) continue; map nResB = nRes; int ncountB = ncount; findChildren(p2,nResB,ncountB); if(ncountB!=0) continue; bool matched2 = true; for(auto const & val : nResB) { if(val.second!=0) { matched2 = false; break; } } if(matched2) { - if(p2.pdgId()==223) - _c_omega->fill(event.weight()); - else if(p2.pdgId()==333) - _c_phi->fill(event.weight()); + if(p2.pid()==223) + _c_omega->fill(); + else if(p2.pid()==333) + _c_phi->fill(); foundOmegaPhi=true; } } } if(found && !foundOmegaPhi) - _c_rho->fill(weight); + _c_rho->fill(); } /// Normalise histograms etc., after the run void finalize() { double fact = crossSection()/nanobarn/sumOfWeights(); for(unsigned int ix=1;ix<5;++ix) { double sigma(0.),error(0.); if(ix==1) { sigma = _c_all->val()*fact; error = _c_all->err()*fact; } else if(ix==2) { sigma = _c_omega->val()*fact; error = _c_omega->err()*fact; } else if(ix==3) { sigma = _c_phi->val()*fact; error = _c_phi->err()*fact; } else if(ix==4) { sigma = _c_rho->val()*fact; error = _c_rho->err()*fact; } Scatter2D temphisto(refData(1, 1, ix)); - Scatter2DPtr mult = bookScatter2D(1, 1, ix); + Scatter2DPtr mult; + book(mult, 1, 1, ix); for (size_t b = 0; b < temphisto.numPoints(); b++) { const double x = temphisto.point(b).x(); pair ex = temphisto.point(b).xErrs(); pair ex2 = ex; if(ex2.first ==0.) ex2. first=0.0001; if(ex2.second==0.) ex2.second=0.0001; if (inRange(sqrtS()/GeV, x-ex2.first, x+ex2.second)) { mult->addPoint(x, sigma, ex, make_pair(error,error)); } else { mult->addPoint(x, 0., ex, make_pair(0.,.0)); } } } } //@} /// @name Histograms //@{ CounterPtr _c_all, _c_omega, _c_phi,_c_rho; //@} }; // The hook for the plugin system DECLARE_RIVET_PLUGIN(SND_2019_I1726419); }