Loading [MathJax]/extensions/tex2jax.js
The Gaudi Framework  v31r0 (aeb156f0)
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
ParticlePropertySvc.cpp
Go to the documentation of this file.
1 // ============================================================================
2 // Include files
3 // ============================================================================
4 // STD&STL
5 // ============================================================================
6 #include "boost/algorithm/string/classification.hpp"
7 #include "boost/algorithm/string/split.hpp"
8 #include "boost/algorithm/string/trim.hpp"
9 #include <cctype>
10 #include <fstream>
11 namespace ba = boost::algorithm;
12 // ============================================================================
13 // GaudiKernel
14 // ============================================================================
17 #include "GaudiKernel/MsgStream.h"
20 #include "GaudiKernel/System.h"
21 // ============================================================================
22 //#include "GaudiKernel/ToStream.h"
23 // ============================================================================
24 // Local
25 // ============================================================================
26 #include "ParticlePropertySvc.h"
27 // ============================================================================
28 namespace Gaudi {
32  DECLARE_COMPONENT( ParticlePropertySvc )
33  // ============================================================================
45  // ============================================================================
46  ParticlePropertySvc::ParticlePropertySvc( const std::string& name, ISvcLocator* svc ) : base_class( name, svc ) {
48  // Redefine the default name:
49  if ( System::getEnv( "PARAMFILESROOT", m_filename.value() ) ) { m_filename.value() += "/data/ParticleTable.txt"; }
50  }
51  // ============================================================================
53  // ============================================================================
56  if ( sc.isFailure() ) { return sc; }
57 
58  sc = setProperties();
59  if ( sc.isFailure() ) {
60  error() << " Could not set the properties " << endmsg;
61  return sc;
62  }
63 
64  m_fileAccess = service( "VFSSvc" );
65  if ( !m_fileAccess ) {
66  error() << " Cannot retrieve the VFS service " << endmsg;
67  return StatusCode::FAILURE;
68  }
69 
70  sc = parse();
71  if ( sc.isFailure() ) {
72  error() << " Could not parse the file " << endmsg;
73  return sc;
74  }
75  if ( !m_particles.empty() ) {
76  sc = addParticles();
77  if ( sc.isFailure() ) {
78  error() << " Could not treat particles! " << endmsg;
79  return sc;
80  }
81  }
82  debug() << "ParticleProperties parsed successfully" << endmsg;
83 
84  debug() << "Access properties" << endmsg;
85  // For debugging purposes print out the size of the internal maps
86  // particle name as key: all particles in database are present here
87  debug() << "NameMap size =" << m_namemap.size() << endmsg;
88  // Geant3 ID as key: all particles in database are present here
89  debug() << "GeantID Map size =" << m_idmap.size() << endmsg;
90  // StdHep ID as key: some particles have no valid StdHep ID
91  debug() << "StdHepID Map size =" << m_stdhepidmap.size() << endmsg;
92  // Pythia ID as key: some particles are not defined in Pythia
93  debug() << "PythiaID Map size =" << m_pythiaidmap.size() << endmsg;
94 
95  if ( !m_replaced.empty() ) {
96  info() << "Properties have been redefined for " << m_replaced.size()
97  << " particles : " << Gaudi::Utils::toString( m_replaced ) << endmsg;
98  }
99 
100  return StatusCode::SUCCESS;
101  }
102  // =============================================================================
104  // =============================================================================
106  if ( !m_other.empty() ) {
107  info() << "Additional Properties have been read from files: " << Gaudi::Utils::toString( m_other ) << endmsg;
108  }
109 
110  if ( !m_replaced.empty() ) {
111  always() << "Properties have been redefined for " << m_replaced.size()
112  << " particles : " << Gaudi::Utils::toString( m_replaced ) << endmsg;
113  }
114 
116 
118  return Service::finalize();
119  }
120  // =============================================================================
122  // =============================================================================
123  StatusCode ParticlePropertySvc::push_back( const std::string& particle, int geantId, int jetsetId, double charge,
124  double mass, double tlife, const std::string& evtName, int pythiaId,
125  double maxWidth ) {
126  //
127  auto i = m_owned.insert( std::make_unique<ParticleProperty>( particle, geantId, jetsetId, charge, mass, tlife,
128  evtName, pythiaId, maxWidth ) );
129  //
130  return i.second ? push_back( i.first->get() ) : StatusCode::FAILURE;
131  }
132  // =============================================================================
134  // =============================================================================
136  if ( !pp ) { return StatusCode::FAILURE; }
137  //
138  { // try to add into Geant(3)ID map
139  const int ID = pp->geantID();
140  // is this already in the map?
141  auto ifind = m_idmap.find( ID );
142  if ( m_idmap.end() != ifind && 0 != m_idmap[ID] ) {
143  diff( ifind->second, pp );
144  m_replaced.insert( m_idmap[ID]->particle() );
145  }
146  // put it into the map
147  m_idmap[ID] = pp;
148  }
149  //
150  { // try to add into Name map
151  const std::string& particle = pp->particle();
152  // is this already in the map?
153  auto ifind = m_namemap.find( particle );
154  if ( ifind != m_namemap.end() && ifind->second ) {
155  diff( ifind->second, pp );
156  m_replaced.insert( ifind->second->particle() );
157  }
158  // put it into the map
159  m_namemap[particle] = pp;
160  }
161  //
162  // add to StdHep map only if StdHep ID different from zero and if
163  // not Cerenkov (StdHep ID = gamma)
164  if ( 0 != pp->jetsetID() && "Tcherenkov" != pp->particle() ) { // try to add into StdHepID map
165  const int ID = pp->jetsetID();
166  // is this already in the map?
167  auto ifind = m_stdhepidmap.find( ID );
168  if ( m_stdhepidmap.end() != ifind && ifind->second ) {
169  diff( ifind->second, pp );
170  m_replaced.insert( ifind->second->particle() );
171  }
172  // put it into the map
173  m_stdhepidmap[ID] = pp;
174  }
175  //
176  // add to Pythia map only if Pythia ID is different from
177  // zero ( StdHep id is always different from zero in this case )
178  if ( 0 != pp->pythiaID() && 0 != pp->jetsetID() && "Tcherenkov" != pp->particle() ) { // try to add into PythiaID
179  // map
180  const int ID = pp->pythiaID();
181  // is this already in the map?
182  auto ifind = m_pythiaidmap.find( ID );
183  if ( m_pythiaidmap.end() != ifind && ifind->second ) {
184  diff( ifind->second, pp );
185  m_replaced.insert( ifind->second->particle() );
186  }
187  // put it into the map
188  m_pythiaidmap[ID] = pp;
189  }
190  //
191  return rebuild();
192  }
193  // =============================================================================
195  // =============================================================================
196  namespace {
197  template <class MAP>
198  void _remove_( MAP& m, const ParticleProperty* pp ) {
199  auto i = std::find_if( m.begin(), m.end(), [&]( typename MAP::const_reference i ) { return i.second == pp; } );
200  if ( i != m.end() ) { m.erase( i ); }
201  }
202  } // namespace
203  // ============================================================================
205  if ( !pp ) { return StatusCode::FAILURE; }
206 
207  _remove_( m_idmap, pp );
208  _remove_( m_namemap, pp );
209  _remove_( m_stdhepidmap, pp );
210  _remove_( m_pythiaidmap, pp );
211  //
212  return rebuild();
213  }
214  // ============================================================================
216  // ============================================================================
218 
219  // parse "the main" file
220  StatusCode sc = parse( m_filename );
221  if ( sc.isFailure() ) { return sc; }
222 
223  // parse "other" files
224  for ( auto& file : m_other ) {
225  sc = parse( file );
226  if ( sc.isFailure() ) { return sc; }
227  }
228 
229  // Now check that the file format was consistent with what parser
230  // expected
231  if ( m_namemap.empty() ) {
232  error() << "Format of input file inconsistent with what expected"
233  << " - Check you are using ParticleData.txt" << endmsg;
234  return StatusCode::FAILURE;
235  }
236 
237  return sc;
238  }
239  // ============================================================================
241  auto infile = ( m_fileAccess ? m_fileAccess->open( file ) : nullptr );
242  if ( !infile ) {
243  error() << "Unable to open properties file : " << file << endmsg;
244  return StatusCode::FAILURE;
245  }
246 
248  info() << "Opened particle properties file : " << file << endmsg;
249 
251  tokens.reserve( 9 );
253  while ( std::getline( *infile, line ) ) {
254  // parse each line of the file (comment lines begin with # in the cdf
255  // file,
256  if ( line.front() == '#' ) continue;
257 
258  tokens.clear();
259  ba::trim_left_if( line, ba::is_space() );
260  ba::split( tokens, line, ba::is_space(), boost::token_compress_on );
261  if ( tokens.size() != 9 ) continue;
262 
263  auto gid = std::stoi( tokens[1] );
264  auto jid = std::stoi( tokens[2] );
265  // Change the particles that do not correspond to a pdg number
266  if ( jid == 0 ) jid = 10000000 * gid;
267 
268  // add a particle property
269  sc = push_back( tokens[0], gid, jid, std::stod( tokens[3] ), std::stod( tokens[4] ) * Gaudi::Units::GeV,
270  std::stod( tokens[5] ) * Gaudi::Units::s, tokens[6], std::stoi( tokens[7] ),
271  std::stod( tokens[8] ) * Gaudi::Units::GeV );
272 
273  if ( sc.isFailure() ) {
274  error() << "Error from ParticlePropertySvc::push_back for particle='" << tokens[0] << "'" << endmsg;
275  }
276  }
277  return StatusCode::SUCCESS;
278  }
279  // ============================================================================
285  // ============================================================================
287  if ( !pp ) { return nullptr; }
288  const int ID = pp->pdgID();
289  const int antiID = -1 * ID;
290  for ( const auto& ap : m_vectpp ) {
291  if ( ap && antiID == ap->pdgID() ) { return ap; } // RETURN
292  };
293  //
294  return pp; // RETURN
295  }
296  // ============================================================================
301  // ============================================================================
303  // initialize particle<-->antiParticle relations
304  for ( auto& pp : m_vectpp ) {
305  if ( !pp ) { continue; } // CONTINUE
306  const ParticleProperty* ap = anti( pp );
307  if ( ap ) { pp->setAntiParticle( ap ); }
308  }
309  return StatusCode::SUCCESS;
310  }
311  // ============================================================================
313  // ============================================================================
314  namespace {
316  template <typename MAP, typename SET>
317  void _load_( MAP& m, SET& result ) {
318  for ( auto i = m.begin(); m.end() != i; ++i ) { result.insert( i->second ); }
319  }
320  } // namespace
321  // ============================================================================
323  std::set<mapped_type> local;
324  m_vectpp.clear();
325  m_vectpp.reserve( m_idmap.size() + 100 );
326  // load information from maps into the set
327  _load_( m_idmap, local );
328  _load_( m_namemap, local );
329  _load_( m_stdhepidmap, local );
330  _load_( m_pythiaidmap, local );
331  // load information from set to the linear container vector
332  std::copy( std::begin( local ), std::end( local ), std::back_inserter( m_vectpp ) );
333  return setAntiParticles();
334  }
335  // ============================================================================
336  // treat additional particles
337  // ============================================================================
339  // loop over all "explicit" particles
340  for ( const auto& item : m_particles ) {
341  std::istringstream input( item );
342  // get the name
343  std::string p_name;
344  int p_geant;
345  int p_jetset;
346  double p_charge;
347  double p_mass;
348  double p_ltime;
349  std::string p_evtgen;
350  int p_pythia;
351  double p_maxwid;
352  if ( input >> p_name >> p_geant >> p_jetset >> p_charge >> p_mass >> p_ltime >> p_evtgen >> p_pythia >>
353  p_maxwid ) {
354  always() << " Add/Modify the particle: "
355  << " name='" << p_name << "'"
356  << " geant=" << p_geant << " jetset=" << p_jetset << " charge=" << p_charge << " mass=" << p_mass
357  << " ltime=" << p_ltime << " evtgen='" << p_evtgen << "'"
358  << " pythia=" << p_pythia << " maxwid=" << p_maxwid << endmsg;
359  //
360  StatusCode sc = push_back( p_name, p_geant, p_jetset, p_charge, p_mass * Gaudi::Units::GeV,
361  p_ltime * Gaudi::Units::s, p_evtgen, p_pythia, p_maxwid * Gaudi::Units::GeV );
362  if ( sc.isFailure() ) { return sc; } // RETURN
363  } else {
364  error() << " could not parse '" << item << "'" << endmsg;
365  return StatusCode::FAILURE; // RETURN
366  }
367  }
368  //
369  return StatusCode::SUCCESS;
370  }
371 // ============================================================================
372 #ifdef __ICC
373 // disable icc remark #1572: floating-point equality and inequality comparisons are unreliable
374 // The comparison are meant
375 # pragma warning( push )
376 # pragma warning( disable : 1572 )
377 #endif
379  //
380  if ( o == n ) { return false; }
381  //
382  auto& log = msgStream();
383  log << l;
384  if ( !o || !n ) {
385  log << MSG::WARNING << " ParticleProperty* point to NULL" << endmsg;
386  return true; // RETURN
387  }
388  //
389  bool result = false;
390  if ( o->particle() != n->particle() ) {
391  result = true;
392  log << " Name:'" << o->particle() << "'/'" << n->particle() << "'";
393  }
394  if ( o->geantID() != n->geantID() ) {
395  result = true;
396  log << " G3ID:" << o->geantID() << "/" << n->geantID() << "'";
397  }
398  if ( o->pdgID() != n->pdgID() ) {
399  result = true;
400  log << " PDGID:" << o->pdgID() << "/" << n->pdgID() << "'";
401  }
402  if ( o->pythiaID() != n->pythiaID() ) {
403  result = true;
404  log << " PYID:" << o->pythiaID() << "/" << n->pythiaID() << "'";
405  }
406  if ( o->charge() != n->charge() ) {
407  result = true;
408  log << " Q:" << o->charge() << "/" << n->charge() << "'";
409  }
410  if ( o->mass() != n->mass() ) {
411  result = true;
412  log << " M:" << o->mass() << "/" << n->mass() << "'";
413  }
414  if ( o->lifetime() != n->lifetime() ) {
415  result = true;
416  log << " T:" << o->lifetime() << "/" << n->lifetime() << "'";
417  }
418  if ( o->evtGenName() != n->evtGenName() ) {
419  result = true;
420  log << " EvtGen:" << o->evtGenName() << "/" << n->evtGenName() << "'";
421  }
422  if ( o->maxWidth() != n->maxWidth() ) {
423  result = true;
424  log << " WMAX:" << o->maxWidth() << "/" << n->maxWidth() << "'";
425  }
426  if ( result ) { log << endmsg; }
427  //
428  return result;
429  }
430 } // namespace Gaudi
431 #ifdef __ICC
432 // re-enable icc remark #1572
433 # pragma warning( pop )
434 #endif
435 // ============================================================================
436 // The END
437 // ============================================================================
int pdgID() const
Get the PDG (= JETSET) ID.
GAUDI_API std::string getEnv(const char *var)
get a particular environment variable (returning "UNKNOWN" if not set)
Definition: System.cpp:380
StatusCode rebuild()
rebuild "the linear container" from the map
StatusCode initialize() override
Definition: Service.cpp:60
T empty(T...args)
T copy(T...args)
The ISvcLocator is the interface implemented by the Service Factory in the Application Manager to loc...
Definition: ISvcLocator.h:25
constexpr double s
int geantID() const
Get the GEANT3 ID.
StatusCode finalize() override
Definition: Service.cpp:164
This service provides access to particle properties.
MsgStream & info() const
shortcut for the method msgStream(MSG::INFO)
double maxWidth() const
Get the max width deviation.
T front(T...args)
A trivial class to hold information about a single particle properties.
std::string toString(const TYPE &obj)
the generic implementation of the type conversion to the string
Definition: ToStream.h:334
StatusCode erase(int geantId) override
Erase a property by geant3 id.
T stod(T...args)
T getline(T...args)
StatusCode finalize() override
Finalise the service.
constexpr static const auto SUCCESS
Definition: StatusCode.h:85
StatusCode initialize() override
Initialise the service.
STL namespace.
MapName m_namemap
Map for particle names.
StatusCode parse()
Parses the file and fill all the maps.
T end(T...args)
SmartIF< IFileAccess > m_fileAccess
bool isFailure() const
Definition: StatusCode.h:130
StatusCode push_back(const std::string &particle, int geantId, int jetsetId, double charge, double mass, double tlife, const std::string &evtName, int pythiaId, double maxWidth) override
Create a new particle property.
STL class.
#define DECLARE_COMPONENT(type)
MsgStream & error() const
shortcut for the method msgStream(MSG::ERROR)
std::set< std::string > m_replaced
int jetsetID() const
Get the JETSET(StdHep) ID.
This class is used for returning status codes from appropriate routines.
Definition: StatusCode.h:50
double lifetime() const
Get the particle lifetime.
constexpr double m
Definition: SystemOfUnits.h:92
T clear(T...args)
StatusCode setProperties()
Method for setting declared properties to the values specified for the job.
Definition: Service.cpp:280
MapStdHepID m_stdhepidmap
Map for StdHep Ids.
dictionary l
Definition: gaudirun.py:517
T insert(T...args)
const ParticleProperty * anti(const ParticleProperty *pp) const
helper (protected) function to find an antiparticle for the given particle ID (StdHepID) ...
T find(T...args)
T size(T...args)
STL class.
MsgStream & debug() const
shortcut for the method msgStream(MSG::DEBUG)
T begin(T...args)
bool diff(const ParticleProperty *o, const ParticleProperty *n, const MSG::Level l=MSG::DEBUG) const
T back_inserter(T...args)
#define SET(x)
Gaudi::Property< std::vector< std::string > > m_other
double charge() const
Get the particle charge.
MsgStream & msgStream() const
Return an uninitialized MsgStream.
constexpr static const auto FAILURE
Definition: StatusCode.h:86
std::set< std::unique_ptr< ParticleProperty > > m_owned
MapID m_idmap
Map for geant IDs.
Gaudi::Property< std::string > m_filename
Gaudi::Property< std::vector< std::string > > m_particles
StatusCode service(const std::string &name, const T *&psvc, bool createIf=true) const
Access a service by name, creating it if it doesn&#39;t already exist.
Definition: Service.h:83
constexpr double GeV
void reset(TYPE *ptr=nullptr)
Set the internal pointer to the passed one disposing of the old one.
Definition: SmartIF.h:86
StatusCode setAntiParticles()
helper (protected) function to set the valid particle<–>antiparticle relations
VectPP m_vectpp
Vector of all particle properties.
const std::string & particle() const
Get the particle name.
MsgStream & always() const
shortcut for the method msgStream(MSG::ALWAYS)
double mass() const
Get the particle mass.
int pythiaID() const
Get the Pythia ID.
T stoi(T...args)
virtual std::unique_ptr< std::istream > open(const std::string &url)=0
Find the URL and returns a unique_ptr to an input stream interface of an object that can be used to r...
Helper functions to set/get the application return code.
Definition: __init__.py:1
MsgStream & endmsg(MsgStream &s)
MsgStream Modifier: endmsg. Calls the output method of the MsgStream.
Definition: MsgStream.h:192
MapPythiaID m_pythiaidmap
Map for Pythia Ids.
T reserve(T...args)
const std::string & evtGenName() const
Get the EvtGen name.