ROOTPWA
particleDataTable.cc
Go to the documentation of this file.
1 
2 //
3 // Copyright 2010
4 //
5 // This file is part of rootpwa
6 //
7 // rootpwa is free software: you can redistribute it and/or modify
8 // it under the terms of the GNU General Public License as published by
9 // the Free Software Foundation, either version 3 of the License, or
10 // (at your option) any later version.
11 //
12 // rootpwa is distributed in the hope that it will be useful,
13 // but WITHOUT ANY WARRANTY; without even the implied warranty of
14 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 // GNU General Public License for more details.
16 //
17 // You should have received a copy of the GNU General Public License
18 // along with rootpwa. If not, see <http://www.gnu.org/licenses/>.
19 //
21 //-------------------------------------------------------------------------
22 // File and Version Information:
23 // $Rev:: $: revision of last commit
24 // $Author:: $: author of last commit
25 // $Date:: $: date of last commit
26 //
27 // Description:
28 // singleton class that manages all particle data
29 //
30 //
31 // Author List:
32 // Boris Grube TUM (original author)
33 //
34 //
35 //-------------------------------------------------------------------------
36 
37 
38 #include <fstream>
39 #include <iomanip>
40 #include <iterator>
41 
42 #include <boost/assign/list_inserter.hpp>
43 
44 #include "libConfigUtils.hpp"
45 #include "reportingUtils.hpp"
46 #include "particleDataTable.h"
47 #include "libconfig.h++"
48 
49 
50 using namespace std;
51 using namespace rpwa;
52 using namespace libconfig;
53 using namespace boost::bimaps;
54 
55 
56 namespace {
57 
58  typedef bimap<string, unsigned int> nameGeantIdBimap;
59  nameGeantIdBimap initNameGeantIdTranslator()
60  {
61  nameGeantIdBimap translator;
62  boost::assign::insert(translator)
63  (nameGeantIdBimap::value_type("gamma", 1))
64  (nameGeantIdBimap::value_type("e+", 2))
65  (nameGeantIdBimap::value_type("e-", 3))
66  (nameGeantIdBimap::value_type("mu+", 5))
67  (nameGeantIdBimap::value_type("mu-", 6))
68  (nameGeantIdBimap::value_type("pi0", 7))
69  (nameGeantIdBimap::value_type("pi+", 8))
70  (nameGeantIdBimap::value_type("pi-", 9))
71  (nameGeantIdBimap::value_type("K_L", 10))
72  (nameGeantIdBimap::value_type("K+", 11))
73  (nameGeantIdBimap::value_type("K-", 12))
74  (nameGeantIdBimap::value_type("n", 13))
75  (nameGeantIdBimap::value_type("p+", 14))
76  (nameGeantIdBimap::value_type("pbar-", 15))
77  (nameGeantIdBimap::value_type("K_S", 16))
78  (nameGeantIdBimap::value_type("eta", 17))
79  (nameGeantIdBimap::value_type("Lambda", 18))
80  (nameGeantIdBimap::value_type("nbar", 25))
81  (nameGeantIdBimap::value_type("Lambdabar", 26))
82  (nameGeantIdBimap::value_type("rho(770)0", 57))
83  (nameGeantIdBimap::value_type("rho(770)+", 58))
84  (nameGeantIdBimap::value_type("rho(770)-", 59))
85  (nameGeantIdBimap::value_type("omega(782)", 60))
86  (nameGeantIdBimap::value_type("eta'(982)", 61))
87  (nameGeantIdBimap::value_type("phi(1020)", 62));
88  return translator;
89  }
90 
91 }
92 
93 
94 particleDataTable particleDataTable::_instance;
95 map<string, particleProperties> particleDataTable::_dataTable;
96 bool particleDataTable::_debug = false;
97 nameGeantIdBimap particleDataTable::_nameGeantIdMap = initNameGeantIdTranslator();
98 
99 
100 string particleDataTable::particleNameFromGeantId(const int id) {
101  nameGeantIdBimap::right_const_iterator it = _nameGeantIdMap.right.find(id);
102  if (it == _nameGeantIdMap.right.end()) {
103  printErr << id << " is unknown GEANT particle ID. returning particle name 'unknown'." << endl;
104  return "unknown";
105  }
106  // make sure charge is correctly put into particle name
107  int charge;
108  const string bareName = particleProperties::chargeFromName(it->second, charge);
109  return particleProperties::nameWithCharge(bareName, charge);
110 }
111 
112 
113 void particleDataTable::geantIdAndChargeFromParticleName(const string& name,
114  int& id,
115  int& charge)
116 {
117  id = 0;
118  const string bareName = particleProperties::chargeFromName(name, charge);
119  nameGeantIdBimap::left_const_iterator it = _nameGeantIdMap.left.find(name);
120  if(it == _nameGeantIdMap.left.end()) {
121  // try again with charge stripped from name
122  it = _nameGeantIdMap.left.find(bareName);
123  }
124  if(it == _nameGeantIdMap.left.end()) {
125  printErr << "particle '" << name << "' is unknown. returning GEANT particle ID 0." << endl;
126  return;
127  }
128  id = it->second;
129 }
130 
131 
132 unsigned int particleDataTable::geantIdFromParticleName(const std::string& name) {
133  int retval;
134  int tmp;
135  geantIdAndChargeFromParticleName(name, retval, tmp);
136  return retval;
137 }
138 
139 
140 bool
141 particleDataTable::isInTable(const string& partName)
142 {
143  if (_dataTable.find(partName) == _dataTable.end()) {
144  return false;
145  } else
146  return true;
147 }
148 
149 
150 const particleProperties*
151 particleDataTable::entry(const string& partName,
152  const bool warnIfNotExistent)
153 {
154  iterator i = _dataTable.find(partName);
155  if (i == _dataTable.end()) {
156  if (warnIfNotExistent)
157  printWarn << "could not find entry for particle '" << partName << "'" << endl;
158  return 0;
159  } else
160  return &(i->second);
161 }
162 
163 
164 vector<const particleProperties*>
165 particleDataTable::entriesMatching(const particleProperties& prototype,
166  const string& sel,
167  const double minMass,
168  const double minMassWidthFactor,
169  const vector<string>& whiteList,
170  const vector<string>& blackList,
171  const particleProperties::decayMode& decay,
172  const bool& forceDecayCheck)
173 {
174  const pair<particleProperties, string> selector(prototype, sel);
175  vector<const particleProperties*> matchingEntries;
176  for (iterator i = _dataTable.begin(); i != _dataTable.end(); ++i) {
177  if (i->second != selector)
178  continue;
179  const particleProperties& partProp = i->second;
180  // limit isobar mass, if minMass > 0
181  if ((minMass > 0) and (partProp.mass() + minMassWidthFactor * partProp.width() < minMass))
182  {
183  if(_debug)printDebug << partProp.name() << " not in mass window " << flush;
184  continue;
185  }
186  // apply white list
187  bool whiteListMatch = (whiteList.size() == 0) ? true : false;
188  for (size_t j = 0; j < whiteList.size(); ++j)
189  if ((partProp.name() == whiteList[j]) or (partProp.bareName() == whiteList[j])) {
190  whiteListMatch = true;
191  break;
192  }
193  if (not whiteListMatch){
194  if(_debug)printDebug << partProp.name() << " not in whitelist " << endl;
195  continue;
196  }
197  // apply black list
198  bool blackListMatch = false;
199  for (size_t j = 0; j < blackList.size(); ++j)
200  if ((partProp.name() == blackList[j]) or (partProp.bareName() == blackList[j])) {
201  blackListMatch = true;
202  break;
203  }
204  if (blackListMatch) {
205  if(_debug)printDebug << partProp.name() << " on blacklist " << endl;
206  continue;
207  }
208  // apply list of decays
209  bool decayMatch = true;
210  if ((decay._daughters.size() > 0) and partProp.nmbDecays() > 0) {
211  if (not partProp.hasDecay(decay))
212  decayMatch = false;
213  } else if (forceDecayCheck)
214  decayMatch = false;
215  if (not decayMatch) {
216  if (_debug)
217  printDebug << partProp.name() << " does not have a decay into " << decay << endl;
218  continue;
219  }
220 
221  if (_debug) {
222  printDebug << "found entry " << partProp.name() << " matching " << prototype
223  << " and '" << sel << "'" << flush;
224  if (minMass > 0)
225  cout << " with mass > " << minMass - minMassWidthFactor * partProp.width() << " GeV";
226  if (whiteList.size() > 0)
227  cout << " ; in white list";
228  if (blackList.size() > 0)
229  cout << " ; not in black list";
230  if (decayMatch)
231  cout << " ; with allowed decay into " << decay << endl;
232  }
233  matchingEntries.push_back(&partProp);
234  }
235  return matchingEntries;
236 }
237 
238 
239 bool
240 particleDataTable::addEntry(const particleProperties& partProp)
241 {
242  const string name = partProp.name();
243  iterator i = _dataTable.find(name);
244  if (i != _dataTable.end()) {
245  printWarn << "trying to add entry for particle '" << name << "' "
246  << "which already exists in table" << endl
247  << " existing entry: " << i->second << endl
248  << " conflicts with: " << partProp << endl
249  << " entry was not added to table." << endl;
250  return false;
251  } else {
252  _dataTable[name] = partProp;
253  if (_debug)
254  printDebug << "added entry for '" << name << "' into particle data table" << endl;
255  return true;
256  }
257 }
258 
259 
260 ostream&
261 particleDataTable::print(ostream& out)
262 {
263  unsigned int countEntries = 0;
264  for (iterator i = begin(); i != end(); ++i) {
265  ++countEntries;
266  out << "entry " << setw(3) << countEntries << ": " << i->second << endl;
267  }
268  return out;
269 }
270 
271 
272 ostream&
273 particleDataTable::dump(ostream& out)
274 {
275  for (iterator i = begin(); i != end(); ++i) {
276  i->second.dump(out);
277  out << endl;
278  }
279  return out;
280 }
281 
282 
283 bool
284 particleDataTable::readFile(const string& fileName)
285 {
286  printInfo << "reading particle data from file '" << fileName << "'" << endl;
287  ifstream file(fileName.c_str());
288  if (not file or not file.good()) {
289  printWarn << "cannot open file '" << fileName << "'" << endl;
290  return false;
291  }
292  return read(file);
293 }
294 
295 
296 bool
297 particleDataTable::read(istream& in)
298 {
299  if (not in or not in.good()) {
300  printWarn << "cannot read from input stream" << endl;
301  return false;
302  }
303  if (_debug)
304  printDebug << "data table has " << nmbEntries() << " entries (before reading)" << endl;
305  unsigned int countEntries = 0;
306  while (in.good()) {
307  particleProperties partProp;
308  if (in >> partProp) {
309  if (addEntry(partProp))
310  ++countEntries;
311  if (not partProp.isItsOwnAntiPart() and addEntry(partProp.antiPartProperties()))
312  ++countEntries;
313  }
314  }
315  printSucc << "read " << countEntries << " new entries into particle data table" << endl;
316  if (_debug)
317  cout << " data table has " << nmbEntries() << " entries (after reading)" << endl;
318  return true;
319 }
320 
321 
322 bool
323 particleDataTable::readDecayModeFile(const string& fileName)
324 {
325  printInfo << "reading particle decay modes from file '" << fileName << "'" << endl;
326  // parse decay mode file
327  Config file;
328  if (not parseLibConfigFile(fileName, file, _debug)) {
329  printWarn << "problems reading file with decay modes '" << fileName << "'. "
330  << "cannot set decay modes." << endl;
331  return false;
332  }
333  const Setting& root = file.getRoot();
334  // find list of particle decays
335  const Setting* partDecayList = findLibConfigList(root, "particleDecayList");
336  if (not partDecayList) {
337  printWarn << "cannot find 'particleDecayList' list. cannot set decay modes." << endl;
338  return false;
339  }
340  const unsigned int nmbEntries = partDecayList->getLength();
341  if (_debug)
342  printDebug << "found particle decay list with " << nmbEntries << " entries in file "
343  << "'" << fileName << "'" << endl;
344  // loop over decay list entries and add decays to particleDataTable
345  unsigned int countEntries = 0;
346  for (unsigned int entryIndex = 0; entryIndex < nmbEntries; ++entryIndex) {
347  const Setting& partDecayEntry = (*partDecayList)[entryIndex];
348  string partName;
349  if (not partDecayEntry.lookupValue("name", partName)) {
350  printWarn << "particle decay list entry [" << entryIndex << "] does not hava a 'name' field. "
351  << "ignoring entry." << endl;
352  continue;
353  }
354  if (_debug)
355  printDebug << "setting decays for '" << partName << "':" << endl;
356  // get list of decay modes for this particle
357  const Setting* partDecays = findLibConfigList(partDecayEntry, "decays");
358  if (not partDecays) {
359  printWarn << "cannot find 'decays' list for " << partName
360  << " particle decay (list entry [" << entryIndex << "]). "
361  << "ignoring entry." << endl;
362  continue;
363  }
364  // lookup particle in database
365  map<string, particleProperties>::iterator dataTableEntry = _dataTable.find(partName);
366  if (dataTableEntry == _dataTable.end()) {
367  printWarn << "could not find particle " << partName << " in data table. "
368  << "ignoring entry." << endl;
369  continue;
370  }
371  particleProperties& particleProp = dataTableEntry->second;
372  // loop over decay modes for this particle
373  const unsigned int nmbDecays = partDecays->getLength();
374  for (unsigned int decayIndex = 0; decayIndex < nmbDecays; ++decayIndex) {
375  // get string array of decay daughters
376  const Setting* decayDaughters = findLibConfigArray((*partDecays)[decayIndex], "products");
377  if (not decayDaughters) {
378  printWarn << "cannot find 'products' array for " << partName
379  << " particle decay list entry [" << decayIndex << "]. "
380  << "ignoring entry." << endl;
381  continue;
382  }
383  if (_debug)
384  cout << " -> ";
385  // loop over decay products
386  const unsigned int nmbDaughters = decayDaughters->getLength();
387  multiset<string> daughters;
388  for (unsigned int daughterIndex = 0; daughterIndex < nmbDaughters; ++daughterIndex) {
389  const string daughterName = (*decayDaughters)[daughterIndex];
390  daughters.insert(daughterName);
391  if (_debug)
392  cout << daughterName << ((daughterIndex < nmbDaughters - 1) ? " " : "");
393  }
394  if(_debug)
395  cout << endl;
396  particleProperties::decayMode decay(daughters);
397  (*partDecays)[decayIndex].lookupValue("L", decay._L);
398  (*partDecays)[decayIndex].lookupValue("S", decay._S);
399  particleProp.addDecayMode(decay);
400  ++countEntries;
401  }
402  } // end loop over decay list entries
403 
404  printSucc << "read " << countEntries << " particle decays into particle data table" << endl;
405  return true;
406 }
Generated by   doxygen 1.8.1.2