|
|
[Rivet-svn] r2060 - in trunk: data/refdata src src/Analyses src/Projections
blackhole at projects.hepforge.org
blackhole at projects.hepforge.org
Mon Nov 16 16:31:19 GMT 2009
Author: eike
Date: Mon Nov 16 16:31:18 2009
New Revision: 2060
Log:
* initial BELLE_2006_S6265367 upload
* add sqrt(s) for asymmetric beams
* fix typo(?) in check for consistent beam set up in Run.cc
Added:
trunk/src/Analyses/BELLE_2006_S6265367.cc
Modified:
trunk/data/refdata/Makefile.am
trunk/src/Analyses/Makefile.am
trunk/src/Projections/Beam.cc
trunk/src/Run.cc
Modified: trunk/data/refdata/Makefile.am
==============================================================================
--- trunk/data/refdata/Makefile.am Mon Nov 16 16:17:29 2009 (r2059)
+++ trunk/data/refdata/Makefile.am Mon Nov 16 16:31:18 2009 (r2060)
@@ -81,4 +81,5 @@
E735_1998_S3905616.aida \
SFM_1984_S1178091.aida \
PDG_HADRON_MULTIPLICITIES.aida \
- PDG_HADRON_MULTIPLICITIES_RATIOS.aida
+ PDG_HADRON_MULTIPLICITIES_RATIOS.aida \
+ BELLE_2006_S6265367.aida
Added: trunk/src/Analyses/BELLE_2006_S6265367.cc
==============================================================================
--- /dev/null 00:00:00 1970 (empty, because file is newly added)
+++ trunk/src/Analyses/BELLE_2006_S6265367.cc Mon Nov 16 16:31:18 2009 (r2060)
@@ -0,0 +1,219 @@
+// -*- C++ -*-
+#include "Rivet/Analysis.hh"
+#include "Rivet/RivetAIDA.hh"
+#include "Rivet/Tools/ParticleIdUtils.hh"
+#include "Rivet/Projections/Beam.hh"
+// #include "Rivet/Projections/FinalState.hh"
+#include "Rivet/Projections/UnstableFinalState.hh"
+// #include "Rivet/Projections/InitialQuarks.hh"
+
+namespace Rivet {
+
+
+ /// @brief BELLE charmed mesons and baryons from fragmentation
+ /// @author Eike von Seggern
+ class BELLE_2006_S6265367 : public Analysis {
+ public:
+
+ BELLE_2006_S6265367(): Analysis("BELLE_2006_S6265367")
+ {
+ setBeams(ELECTRON, POSITRON);
+ setNeedsCrossSection(true);
+ }
+
+
+ void analyze(const Event& e) {
+ //
+ // TODO: apply BELLE hadron selection cuts
+ //
+
+ const double weight = e.weight();
+
+ // Loop through unstable FS particles and look for charmed mesons/baryons
+ const UnstableFinalState& ufs = applyProjection<UnstableFinalState>(e, "UFS");
+
+ const Beam beamproj = applyProjection<Beam>(e, "Beams");
+ const ParticlePair& beams = beamproj.beams();
+ const double s = beamproj.sqrtS()*beamproj.sqrtS();
+
+ // TODO: implement sqrtS() for asymm. beams in beam projection
+ const bool onresonance = true;
+ // const bool onresonance = fuzzyEquals(beamproj.sqrtS(), 10.58, 1E-4);
+
+ // Particle masses from PDGlive (online 16. Nov. 2009).
+ foreach (const Particle& p, ufs.particles()) {
+ // TODO: Data is not corrected for branching fractions.
+ // TODO: transform vectors to e+e- rest frame
+
+ double xp = 0.0;
+ double mH2 = 0.0;
+ const double mom = p.momentum().vector3().mod();
+
+ const int PdgId = abs(p.pdgId());
+ getLog() << Log::DEBUG << "pdgID = " << PdgId << " mom = " << mom << endl;
+
+ switch (PdgId) {
+
+ case 421:
+ getLog() << Log::DEBUG << "D0 found" << endl;
+ mH2 = 3.47763; // 1.86484^2
+ xp = mom/sqrt(s/4.0 - mH2);
+ if (onresonance)
+ _histXpD0_R->fill(xp, weight);
+ else
+ _histXpD0_C->fill(xp, weight);
+ break;
+
+ case 411:
+ getLog() << Log::DEBUG << "D+ found" << endl;
+ mH2 = 3.49547; // 1.86962^2
+ xp = mom/sqrt(s/4.0 - mH2);
+ if (onresonance)
+ _histXpDplus_R->fill(xp, weight);
+ else
+ _histXpDplus_C->fill(xp, weight);
+ break;
+
+ case 431:
+ getLog() << Log::DEBUG << "D+_s found" << endl;
+ mH2 = 3.87495; // 1.96849^2
+ xp = mom/sqrt(s/4.0 - mH2);
+ if (onresonance)
+ _histXpDplus_s_R->fill(xp, weight);
+ else
+ _histXpDplus_s_C->fill(xp, weight);
+ break;
+
+ case 4122:
+ getLog() << Log::DEBUG << "Lambda_c found" << endl;
+ mH2 = 5.22780; // 2.28646^2
+ xp = mom/sqrt(s/4.0 - mH2);
+ if (onresonance)
+ _histXpLambda_c_R->fill(xp, weight);
+ else
+ _histXpLambda_c_C->fill(xp, weight);
+ break;
+
+ case 413: {
+ getLog() << Log::DEBUG << "D*+ found" << endl;
+ mH2 = 4.04119; // 2.01027^2
+ xp = mom/sqrt(s/4.0 - mH2);
+
+ const GenVertex* dv = p.genParticle().end_vertex();
+ bool D0decay(false), Pi0decay(false), Piplusdecay(false), Dplusdecay(false);
+
+ for (GenVertex::particles_out_const_iterator
+ pp = dv->particles_out_const_begin();
+ pp != dv->particles_out_const_end(); ++pp) {
+ if (abs((*pp)->pdg_id()) == 421) {
+ D0decay = true;
+ } else if (abs((*pp)->pdg_id()) == 411) {
+ Dplusdecay = true;
+ } else if (abs((*pp)->pdg_id()) == 111) {
+ Pi0decay = true;
+ } else if (abs((*pp)->pdg_id()) == 211) {
+ Piplusdecay = true;
+ }
+ }
+
+ if (D0decay && Piplusdecay) {
+ if (onresonance)
+ _histXpDstarplus2D0_R->fill(xp, weight);
+ else
+ _histXpDstarplus2D0_C->fill(xp, weight);
+ } else if (Dplusdecay && Pi0decay) {
+ if (onresonance)
+ _histXpDstarplus2Dplus_R->fill(xp, weight);
+ else
+ _histXpDstarplus2Dplus_C->fill(xp, weight);
+ } else {
+ getLog() << Log::WARN << "Unexpected D* decay!" << endl;
+ }
+ break;
+ }
+
+ case 423:
+ getLog() << Log::DEBUG << "D*0 found" << endl;
+ mH2 = 4.02793; // 2.00697**2
+ xp = mom/sqrt(s/4.0 - mH2);
+ getLog() << Log::DEBUG << "xp = " << xp << endl;
+ if (onresonance)
+ _histXpDstar0_R->fill(xp, weight);
+ else
+ _histXpDstar0_C->fill(xp, weight);
+ break;
+ }
+
+ }
+ } // analyze
+
+
+ void finalize() {
+ normalize(_histXpDstarplus2D0_R, crossSection());
+ normalize(_histXpD0_R, crossSection());
+ normalize(_histXpDplus_R, crossSection());
+ normalize(_histXpDplus_s_R, crossSection());
+ normalize(_histXpLambda_c_R, crossSection());
+ normalize(_histXpDstarplus2Dplus_R, crossSection());
+ normalize(_histXpDstar0_R, crossSection());
+
+ normalize(_histXpDstarplus2D0_C, crossSection());
+ normalize(_histXpD0_C, crossSection());
+ normalize(_histXpDplus_C, crossSection());
+ normalize(_histXpDplus_s_C, crossSection());
+ normalize(_histXpLambda_c_C, crossSection());
+ normalize(_histXpDstarplus2Dplus_C, crossSection());
+ normalize(_histXpDstar0_C, crossSection());
+ } // finalize
+
+
+ void init() {
+ addProjection(Beam(), "Beams");
+ addProjection(UnstableFinalState(-1.3170, 1.9008), "UFS");
+
+ // histograms for continuum data (sqrt(s) = 10.52 GeV)
+ _histXpDstarplus2D0_C = bookHistogram1D(2, 1, 1);
+ _histXpD0_C = bookHistogram1D(3, 1, 1);
+ _histXpDplus_C = bookHistogram1D(4, 1, 1);
+ _histXpDplus_s_C = bookHistogram1D(5, 1, 1);
+ _histXpLambda_c_C = bookHistogram1D(6, 1, 1);
+ _histXpDstarplus2Dplus_C = bookHistogram1D(7, 1, 1);
+ _histXpDstar0_C = bookHistogram1D(8, 1, 1);
+
+ // histograms for on-resonance data (sqrt(s) = 10.58 GeV)
+ _histXpDstarplus2D0_R = bookHistogram1D(9, 1, 1);
+ _histXpD0_R = bookHistogram1D(10, 1, 1);
+ _histXpDplus_R = bookHistogram1D(11, 1, 1);
+ _histXpDplus_s_R = bookHistogram1D(12, 1, 1);
+ _histXpLambda_c_R = bookHistogram1D(13, 1, 1);
+ _histXpDstarplus2Dplus_R = bookHistogram1D(14, 1, 1);
+ _histXpDstar0_R = bookHistogram1D(15, 1, 1);
+
+ } // init
+
+ private:
+
+ //@{
+ /// Histograms
+ // histograms for continuum data (sqrt(s) = 10.52 GeV)
+ AIDA::IHistogram1D* _histXpDstarplus2D0_C;
+ AIDA::IHistogram1D* _histXpD0_C;
+ AIDA::IHistogram1D* _histXpDplus_C;
+ AIDA::IHistogram1D* _histXpDplus_s_C;
+ AIDA::IHistogram1D* _histXpLambda_c_C;
+ AIDA::IHistogram1D* _histXpDstarplus2Dplus_C;
+ AIDA::IHistogram1D* _histXpDstar0_C;
+
+ // histograms for on-resonance data (sqrt(s) = 10.58 GeV)
+ AIDA::IHistogram1D* _histXpDstarplus2D0_R;
+ AIDA::IHistogram1D* _histXpD0_R;
+ AIDA::IHistogram1D* _histXpDplus_R;
+ AIDA::IHistogram1D* _histXpDplus_s_R;
+ AIDA::IHistogram1D* _histXpLambda_c_R;
+ AIDA::IHistogram1D* _histXpDstarplus2Dplus_R;
+ AIDA::IHistogram1D* _histXpDstar0_R;
+ //@}
+ };
+
+ AnalysisBuilder<BELLE_2006_S6265367> plugin_BELLE_2006_S6265367;
+}
Modified: trunk/src/Analyses/Makefile.am
==============================================================================
--- trunk/src/Analyses/Makefile.am Mon Nov 16 16:17:29 2009 (r2059)
+++ trunk/src/Analyses/Makefile.am Mon Nov 16 16:31:18 2009 (r2060)
@@ -105,7 +105,8 @@
PDG_Hadron_Multiplicities_Ratios.cc \
SFM_1984_S1178091.cc \
E735_1998_S3905616.cc \
- JADE_OPAL_2000_S4300807.cc
+ JADE_OPAL_2000_S4300807.cc \
+ BELLE_2006_S6265367.cc
lib_LTLIBRARIES += RivetMCAnalyses.la
Modified: trunk/src/Projections/Beam.cc
==============================================================================
--- trunk/src/Projections/Beam.cc Mon Nov 16 16:17:29 2009 (r2059)
+++ trunk/src/Projections/Beam.cc Mon Nov 16 16:31:18 2009 (r2060)
@@ -30,7 +30,11 @@
sqrts = fabs(mom1) + fabs(mom2);
} else {
/// @todo Implement general sqrt(s) for asymmetric beams... requires particle masses.
- throw Error("Asymmetric beams... calculation of sqrt(S) not yet implemented");
+ getLog() << Log::WARN << "Asymmetric beams: mass treatment still to be implemented!" << endl;
+ const double E1 = beams().first.momentum().E();
+ const double E2 = beams().second.momentum().E();
+ sqrts = (E1+E2)*(E1+E2) - (mom1+mom2)*(mom1+mom2);
+ sqrts = sqrt(sqrts);
}
getLog() << Log::DEBUG << "sqrt(s) = " << sqrts/GeV << " GeV" << endl;
return sqrts;
Modified: trunk/src/Run.cc
==============================================================================
--- trunk/src/Run.cc Mon Nov 16 16:17:29 2009 (r2059)
+++ trunk/src/Run.cc Mon Nov 16 16:31:18 2009 (r2060)
@@ -72,7 +72,7 @@
Log::getLog("Rivet.Run") << Log::INFO << "First event beams: "
<< this->beams() << " @ " << this->sqrtS()/GeV << " GeV" << endl;
} else {
- if (_beams != _beams || !fuzzyEquals(sqrts, sqrtS())) {
+ if (beams != _beams || !fuzzyEquals(sqrts, sqrtS())) {
Log::getLog("Rivet.Run") << Log::ERROR << "Event beams mismatch: "
<< beams << " @ " << sqrts/GeV << " GeV" << " vs. first beams "
<< this->beams() << " @ " << this->sqrtS()/GeV << " GeV" << endl;
More information about the Rivet-svn
mailing list
|