11 #include <generators/touschek/TouschekReaderTURTLE.h>
13 #include <framework/logging/Logger.h>
14 #include <framework/gearbox/Unit.h>
15 #include <mdst/dataobjects/MCParticle.h>
17 #include <boost/foreach.hpp>
18 #include <boost/lexical_cast.hpp>
19 #include <boost/algorithm/string.hpp>
20 #include <boost/tokenizer.hpp>
26 TouschekReaderTURTLE::TouschekReaderTURTLE(TGeoHMatrix* transMatrix,
int pdg): m_transMatrix(transMatrix), m_pdg(pdg),
42 if (!
m_input)
throw(TouschekCouldNotOpenFileError() << filename);
49 B2ERROR(
"The transformation matrix is NULL !");
54 typedef boost::tokenizer<boost::char_separator<char> > tokenizer;
55 boost::char_separator<char> sep(
",; \t");
57 double particlePosTouschek[3] = {0.0, 0.0, 0.0};
58 double particlePosGeant4[3] = {0.0, 0.0, 0.0};
59 double particleMomTouschek[3] = {0.0, 0.0, 0.0};
60 double particleMomGeant4[3] = {0.0, 0.0, 0.0};
66 boost::trim(currLine);
67 if (currLine ==
"")
continue;
69 tokenizer tokens(currLine, sep);
72 BOOST_FOREACH(
const string & tok, tokens) {
75 B2WARNING(
"This Touschek file has more than 7 fields ! Only the first 7 fields were read.");
78 fields[index] = boost::lexical_cast<double>(tok);
80 }
catch (boost::bad_lexical_cast&) {
81 throw (TouschekConvertFieldError() <<
m_lineNum << index << tok);
87 particlePosTouschek[0] = fields[2] *
Unit::m;
88 particlePosTouschek[1] = -fields[3] *
Unit::m;
89 if (
m_pdg == -11) particlePosTouschek[2] = 400.0;
90 else if (
m_pdg == 11) particlePosTouschek[2] = -400.0;
91 else { B2WARNING(
"m_pdg is not 11/-11");
break;}
93 m_transMatrix->LocalToMaster(particlePosTouschek, particlePosGeant4);
97 particleMomTouschek[0] = fields[4];
98 particleMomTouschek[1] = -fields[5];
99 if (
m_pdg == -11) particleMomTouschek[2] = -fields[6];
100 else if (
m_pdg == 11) particleMomTouschek[2] = fields[6];
101 else { B2WARNING(
"m_pdg is not 11/-11");
break;}
103 m_transMatrix->LocalToMasterVect(particleMomTouschek, particleMomGeant4);
105 double totalMom2 = particleMomGeant4[0] * particleMomGeant4[0];
106 totalMom2 += particleMomGeant4[1] * particleMomGeant4[1];
107 totalMom2 += particleMomGeant4[2] * particleMomGeant4[2];
112 particle.setPDG(
m_pdg);
113 particle.setMassFromPDG();
114 particle.setMomentum(TVector3(particleMomGeant4));
115 particle.setProductionVertex(TVector3(particlePosGeant4));
116 particle.setEnergy(sqrt(totalMom2 + particle.getMass()*particle.getMass()));
117 particle.setProductionTime(0.0);
118 particle.setValidVertex(
true);
121 }
while ((!
m_input.eof()) && ((number > -1) && (numParticles < number)));