Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TDatabasePDG.cxx
Go to the documentation of this file.
1// @(#)root/eg:$Id$
2// Author: Pasha Murat 12/02/99
3
4/*************************************************************************
5 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12#include "RConfigure.h"
13#include "TROOT.h"
14#include "TEnv.h"
15#include "THashList.h"
16#include "TExMap.h"
17#include "TSystem.h"
18#include "TDatabasePDG.h"
19#include "TDecayChannel.h"
20#include "TParticlePDG.h"
21#include <cstdlib>
22
23
24/** \class TDatabasePDG
25 \ingroup eg
26
27Particle database manager class
28
29This manager creates a list of particles which by default is
30initialised from with the constants used by PYTHIA6 (plus some
31other particles added). See definition and the format of the default
32particle list in $ROOTSYS/etc/pdg_table.txt
33
34There are 2 ways of redefining the name of the file containing the
35particle properties
36
37-# One can define the name in .rootrc file:
38
39Root.DatabasePDG: $(HOME)/my_pdg_table.txt
40
41-# One can use TDatabasePDG::ReadPDGTable method explicitly:
42
43 - TDatabasePDG *pdg = new TDatabasePDG();
44 - pdg->ReadPDGtable(filename)
45
46See TParticlePDG for the description of a static particle properties.
47See TParticle for the description of a dynamic particle particle.
48
49The current default pdg_table file displays lifetime 0 for some unstable particles.
50
51*/
52
53
54////////////////////////////////////////////////////////////////////////////////
55/// Static function holding the instance.
56
58{
59 static TDatabasePDG* fgInstance = nullptr;
60 return &fgInstance;
61}
62
63////////////////////////////////////////////////////////////////////////////////
64/// Create PDG database. Initialization of the DB has to be done via explicit
65/// call to ReadDataBasePDG (also done by GetParticle methods)
66
67TDatabasePDG::TDatabasePDG(): TNamed("PDGDB","The PDG particle data base")
68{
69 fParticleList = 0;
70 fPdgMap = 0;
72 auto fgInstance = GetInstancePtr();
73 if (*fgInstance != nullptr) {
74 Warning("TDatabasePDG", "object already instantiated");
75 } else {
76 *fgInstance = this;
77 gROOT->GetListOfSpecials()->Add(this);
78 }
79}
80
81////////////////////////////////////////////////////////////////////////////////
82/// Cleanup the PDG database.
83
85{
86 if (fParticleList) {
88 delete fParticleList; // this deletes all objects in the list
89 if (fPdgMap) delete fPdgMap;
90 }
91 // classes do not own particles...
92 if (fListOfClasses) {
94 delete fListOfClasses;
95 }
96 if (gROOT && !gROOT->TestBit(TObject::kInvalidObject))
97 gROOT->GetListOfSpecials()->Remove(this);
98 auto fgInstance = GetInstancePtr();
99 *fgInstance = nullptr;
100}
101
102////////////////////////////////////////////////////////////////////////////////
103///static function
104
106{
107 auto fgInstance = GetInstancePtr();
109 if (*fgInstance == nullptr) {
111 if (!*fgInstance) {
112 // Constructor creates a new instance, inits fgInstance.
113 new TDatabasePDG();
114 }
115 }
116 return *fgInstance;
117}
118
119////////////////////////////////////////////////////////////////////////////////
120/// Build fPdgMap mapping pdg-code to particle.
121///
122/// Initial size is set so as to be able to hold at least 600
123/// particles: 521 in default table, ALICE adds 54 more.
124/// To be revisited after LHC discovers SUSY.
125
127{
128 // It is preferrable to waste some work in presence of high contention
129 // for the benefit of reducing the critical section to the bare minimum
130 auto pdgMap = new TExMap(4 * std::max(600, fParticleList->GetEntries()) / 3 + 3);
131 TIter next(fParticleList);
133 while ((p = (TParticlePDG*)next())) {
134 pdgMap->Add((Long64_t)p->PdgCode(), (Long64_t)p);
135 }
136
138 if (!fPdgMap) {
139 fPdgMap = pdgMap;
140 } else {
141 delete pdgMap;
142 }
143}
144
145////////////////////////////////////////////////////////////////////////////////
146///
147/// Particle definition normal constructor. If the particle is set to be
148/// stable, the decay width parameter does have no meaning and can be set to
149/// any value. The parameters granularity, LowerCutOff and HighCutOff are
150/// used for the construction of the mean free path look up tables. The
151/// granularity will be the number of logwise energy points for which the
152/// mean free path will be calculated.
153///
154
155TParticlePDG* TDatabasePDG::AddParticle(const char *name, const char *title,
156 Double_t mass, Bool_t stable,
158 const char* ParticleClass,
160 Int_t Anti,
161 Int_t TrackingCode)
162{
164
165 if (old) {
166 Warning("AddParticle", "Particle with PDGcode=%d already defined", PDGcode);
167 return 0;
168 }
169
170 // Fix ROOT-6889. If the antiparticle is stable,
171 if (PDGcode < 0 && stable) {
172 // We verify the antiparticle is stable and this is not an artifact of the PDG Table
174 stable = nonAntiParticle->Stable();
175 }
176
177 TParticlePDG* p = new TParticlePDG(name, title, mass, stable, width,
178 charge, ParticleClass, PDGcode, Anti,
179 TrackingCode);
181 if (fPdgMap)
183
184 TParticleClassPDG* pclass = GetParticleClass(ParticleClass);
185
186 if (!pclass) {
187 pclass = new TParticleClassPDG(ParticleClass);
189 }
190
191 pclass->AddParticle(p);
192
193 return p;
194}
195
196////////////////////////////////////////////////////////////////////////////////
197/// assuming particle has already been defined
198
200{
201 TParticlePDG* old = GetParticle(PdgCode);
202
203 if (old) {
204 Warning("AddAntiParticle", "Can't redefine parameters");
205 return NULL;
206 }
207
208 Int_t pdg_code = abs(PdgCode);
210
211 if (!p) {
212 Warning("AddAntiParticle", "Particle with pdg code %d not known", pdg_code);
213 return NULL;
214 }
215
217 Name,
218 p->Mass(),
219 1,
220 p->Width(),
221 -p->Charge(),
222 p->ParticleClass(),
223 PdgCode,
224 1,
225 p->TrackingCode());
226 return ap;
227}
228
229
230////////////////////////////////////////////////////////////////////////////////
231///
232/// Get a pointer to the particle object according to the name given
233///
234
236{
237 {
239 if (fParticleList == 0)
240 ((TDatabasePDG *)this)->ReadPDGTableImpl("", false);
241 }
243// if (!def) {
244// Error("GetParticle","No match for %s exists!",name);
245// }
246 return def;
247}
248
249////////////////////////////////////////////////////////////////////////////////
250///
251/// Get a pointer to the particle object according to the MC code number
252///
253
255{
256 {
258 if (fParticleList == 0)
259 ((TDatabasePDG *)this)->ReadPDGTableImpl("", false);
260 if (fPdgMap == 0)
261 BuildPdgMap();
262 }
263
265}
266
267////////////////////////////////////////////////////////////////////////////////
268/// Print contents of PDG database.
269
271{
272 if (fParticleList == 0) ((TDatabasePDG*)this)->ReadPDGTable();
273
274 TIter next(fParticleList);
276 while ((p = (TParticlePDG *)next())) {
277 p->Print(option);
278 }
279}
280
281////////////////////////////////////////////////////////////////////////////////
282/// Converts Geant3 particle codes to PDG convention. (Geant4 uses
283/// PDG convention already)
284/// Source: BaBar User Guide, Neil I. Geddes,
285///
286/// see <A href="http://www.slac.stanford.edu/BFROOT/www/Computing/Environment/NewUser/htmlbug/node51.html"> Conversion table</A>
287///
288/// with some fixes by PB, marked with (PB) below. Checked against
289/// PDG listings from 2000.
290///
291/// Paul Balm, Nov 19, 2001
292
294{
295
296
297 switch(Geant3number) {
298
299 case 1 : return 22; // photon
300 case 25 : return -2112; // anti-neutron
301 case 2 : return -11; // e+
302 case 26 : return -3122; // anti-Lambda
303 case 3 : return 11; // e-
304 case 27 : return -3222; // Sigma-
305 case 4 : return 12; // e-neutrino (NB: flavour undefined by Geant)
306 case 28 : return -3212; // Sigma0
307 case 5 : return -13; // mu+
308 case 29 : return -3112; // Sigma+ (PB)*/
309 case 6 : return 13; // mu-
310 case 30 : return -3322; // Xi0
311 case 7 : return 111; // pi0
312 case 31 : return -3312; // Xi+
313 case 8 : return 211; // pi+
314 case 32 : return -3334; // Omega+ (PB)
315 case 9 : return -211; // pi-
316 case 33 : return -15; // tau+
317 case 10 : return 130; // K long
318 case 34 : return 15; // tau-
319 case 11 : return 321; // K+
320 case 35 : return 411; // D+
321 case 12 : return -321; // K-
322 case 36 : return -411; // D-
323 case 13 : return 2112; // n
324 case 37 : return 421; // D0
325 case 14 : return 2212; // p
326 case 38 : return -421; // D0
327 case 15 : return -2212; // anti-proton
328 case 39 : return 431; // Ds+
329 case 16 : return 310; // K short
330 case 40 : return -431; // anti Ds-
331 case 17 : return 221; // eta
332 case 41 : return 4122; // Lamba_c+
333 case 18 : return 3122; // Lambda
334 case 42 : return 24; // W+
335 case 19 : return 3222; // Sigma+
336 case 43 : return -24; // W-
337 case 20 : return 3212; // Sigma0
338 case 44 : return 23; // Z
339 case 21 : return 3112; // Sigma-
340 case 45 : return 0; // deuteron
341 case 22 : return 3322; // Xi0
342 case 46 : return 0; // triton
343 case 23 : return 3312; // Xi-
344 case 47 : return 0; // alpha
345 case 24 : return 3334; // Omega- (PB)
346 case 48 : return 0; // G nu ? PDG ID 0 is undefined
347
348 default : return 0;
349
350 }
351}
352
353////////////////////////////////////////////////////////////////////////////////
354/// Converts pdg code to geant3 id
355
357{
358 switch(pdgNumber) {
359
360 case 22 : return 1; // photon
361 case -2112 : return 25; // anti-neutron
362 case -11 : return 2; // e+
363 case -3122 : return 26; // anti-Lambda
364 case 11 : return 3; // e-
365 case -3222 : return 27; // Sigma-
366 case 12 : return 4; // e-neutrino (NB: flavour undefined by Geant)
367 case -3212 : return 28; // Sigma0
368 case -13 : return 5; // mu+
369 case -3112 : return 29; // Sigma+ (PB)*/
370 case 13 : return 6; // mu-
371 case -3322 : return 30; // Xi0
372 case 111 : return 7; // pi0
373 case -3312 : return 31; // Xi+
374 case 211 : return 8; // pi+
375 case -3334 : return 32; // Omega+ (PB)
376 case -211 : return 9; // pi-
377 case -15 : return 33; // tau+
378 case 130 : return 10; // K long
379 case 15 : return 34; // tau-
380 case 321 : return 11; // K+
381 case 411 : return 35; // D+
382 case -321 : return 12; // K-
383 case -411 : return 36; // D-
384 case 2112 : return 13; // n
385 case 421 : return 37; // D0
386 case 2212 : return 14; // p
387 case -421 : return 38; // D0
388 case -2212 : return 15; // anti-proton
389 case 431 : return 39; // Ds+
390 case 310 : return 16; // K short
391 case -431 : return 40; // anti Ds-
392 case 221 : return 17; // eta
393 case 4122 : return 41; // Lamba_c+
394 case 3122 : return 18; // Lambda
395 case 24 : return 42; // W+
396 case 3222 : return 19; // Sigma+
397 case -24 : return 43; // W-
398 case 3212 : return 20; // Sigma0
399 case 23 : return 44; // Z
400 case 3112 : return 21; // Sigma-
401 case 3322 : return 22; // Xi0
402 case 3312 : return 23; // Xi-
403 case 3334 : return 24; // Omega- (PB)
404
405 default : return 0;
406
407 }
408}
409
410////////////////////////////////////////////////////////////////////////////////
411///
412/// Converts the ISAJET Particle number into the PDG MC number
413///
414
416{
417 switch (isaNumber) {
418 case 1 : return 2; // UP .30000E+00 .67
419 case -1 : return -2; // UB .30000E+00 -.67
420 case 2 : return 1; // DN .30000E+00 -.33
421 case -2 : return -1; // DB .30000E+00 .33
422 case 3 : return 3; // ST .50000E+00 -.33
423 case -3 : return -3; // SB .50000E+00 .33
424 case 4 : return 4; // CH .16000E+01 .67
425 case -4 : return -4; // CB .16000E+01 -.67
426 case 5 : return 5; // BT .49000E+01 -.33
427 case -5 : return -5; // BB .49000E+01 .33
428 case 6 : return 6; // TP .17500E+03 .67
429 case -6 : return -6; // TB .17500E+03 -.67
430 case 9 : return 21; // GL 0. 0.00
431 case 80 : return 24; // W+ SIN2W=.23 1.00
432 case -80 : return -24; // W- SIN2W=.23 -1.00
433 case 90 : return 23; // Z0 SIN2W=.23 0.00
434 case 230 : return 311; // K0 .49767E+00 0.00
435 case -230 : return -311; // AK0 .49767E+00 0.00
436 case 330 : return 331; // ETAP .95760E+00 0.00
437 case 340 : return 0; // F- .20300E+01 -1.00
438 case -340 : return 0; // F+ .20300E+01 1.00
439 case 440 : return 441; // ETAC .29760E+01 0.00
440 case 111 : return 113; // RHO0 .77000E+00 0.00
441 case 121 : return 213; // RHO+ .77000E+00 1.00
442 case -121 : return -213; // RHO- .77000E+00 -1.00
443 case 221 : return 223; // OMEG .78260E+00 0.00
444 case 131 : return 323; // K*+ .88810E+00 1.00
445 case -131 : return -323; // K*- .88810E+00 -1.00
446 case 231 : return 313; // K*0 .89220E+00 0.00
447 case -231 : return -313; // AK*0 .89220E+00 0.00
448 case 331 : return 333; // PHI .10196E+01 0.00
449 case -140 : return 421; // D0
450 case 140 : return -421; // D0 bar
451 case 141 : return -423; // AD*0 .20060E+01 0.00
452 case -141 : return 423; // D*0 .20060E+01 0.00
453 case -240 : return -411; // D+
454 case 240 : return 411; // D-
455 case 241 : return -413; // D*- .20086E+01 -1.00
456 case -241 : return 413; // D*+ .20086E+01 1.00
457 case 341 : return 0; // F*- .21400E+01 -1.00
458 case -341 : return 0; // F*+ .21400E+01 1.00
459 case 441 : return 443; // JPSI .30970E+01 0.00
460
461 // B-mesons, Bc still missing
462 case 250 : return 511; // B0
463 case -250 : return -511; // B0 bar
464 case 150 : return 521; // B+
465 case -150 : return -521; // B-
466 case 350 : return 531; // Bs 0
467 case -350 : return -531; // Bs bar
468 case 351 : return 533; // Bs* 0
469 case -351 : return -533; // Bs* bar
470 case 450 : return 541; // Bc +
471 case -450 : return -541; // Bc bar
472
473 case 1140 : return 4222; // SC++ .24300E+01 2.00
474 case -1140 : return -4222; // ASC-- .24300E+01 -2.00
475 case 1240 : return 4212; // SC+ .24300E+01 1.00
476 case -1240 : return -4212; // ASC- .24300E+01 -1.00
477 case 2140 : return 4122; // LC+ .22600E+01 1.00
478 case -2140 : return -4122; // ALC- .22600E+01 -1.00
479 case 2240 : return 4112; // SC0 .24300E+01 0.00
480 case -2240 : return -4112; // ASC0 .24300E+01 0.00
481 case 1340 : return 0; // USC. .25000E+01 1.00
482 case -1340 : return 0; // AUSC. .25000E+01 -1.00
483 case 3140 : return 0; // SUC. .24000E+01 1.00
484 case -3140 : return 0; // ASUC. .24000E+01 -1.00
485 case 2340 : return 0; // DSC. .25000E+01 0.00
486 case -2340 : return 0; // ADSC. .25000E+01 0.00
487 case 3240 : return 0; // SDC. .24000E+01 0.00
488 case -3240 : return 0; // ASDC. .24000E+01 0.00
489 case 3340 : return 0; // SSC. .26000E+01 0.00
490 case -3340 : return 0; // ASSC. .26000E+01 0.00
491 case 1440 : return 0; // UCC. .35500E+01 2.00
492 case -1440 : return 0; // AUCC. .35500E+01 -2.00
493 case 2440 : return 0; // DCC. .35500E+01 1.00
494 case -2440 : return 0; // ADCC. .35500E+01 -1.00
495 case 3440 : return 0; // SCC. .37000E+01 1.00
496 case -3440 : return 0; // ASCC. .37000E+01 -1.00
497 case 1111 : return 2224; // DL++ .12320E+01 2.00
498 case -1111 : return -2224; // ADL-- .12320E+01 -2.00
499 case 1121 : return 2214; // DL+ .12320E+01 1.00
500 case -1121 : return -2214; // ADL- .12320E+01 -1.00
501 case 1221 : return 2114; // DL0 .12320E+01 0.00
502 case -1221 : return -2114; // ADL0 .12320E+01 0.00
503 case 2221 : return 1114; // DL- .12320E+01 -1.00
504 case -2221 : return -1114; // ADL+ .12320E+01 1.00
505 case 1131 : return 3224; // S*+ .13823E+01 1.00
506 case -1131 : return -3224; // AS*- .13823E+01 -1.00
507 case 1231 : return 3214; // S*0 .13820E+01 0.00
508 case -1231 : return -3214; // AS*0 .13820E+01 0.00
509 case 2231 : return 3114; // S*- .13875E+01 -1.00
510 case -2231 : return -3114; // AS*+ .13875E+01 1.00
511 case 1331 : return 3324; // XI*0 .15318E+01 0.00
512 case -1331 : return -3324; // AXI*0 .15318E+01 0.00
513 case 2331 : return 3314; // XI*- .15350E+01 -1.00
514 case -2331 : return -3314; // AXI*+ .15350E+01 1.00
515 case 3331 : return 3334; // OM- .16722E+01 -1.00
516 case -3331 : return -3334; // AOM+ .16722E+01 1.00
517 case 1141 : return 0; // UUC* .26300E+01 2.00
518 case -1141 : return 0; // AUUC* .26300E+01 -2.00
519 case 1241 : return 0; // UDC* .26300E+01 1.00
520 case -1241 : return 0; // AUDC* .26300E+01 -1.00
521 case 2241 : return 0; // DDC* .26300E+01 0.00
522 case -2241 : return 0; // ADDC* .26300E+01 0.00
523 case 1341 : return 0; // USC* .27000E+01 1.00
524 case -1341 : return 0; // AUSC* .27000E+01 -1.00
525 case 2341 : return 0; // DSC* .27000E+01 0.00
526 case -2341 : return 0; // ADSC* .27000E+01 0.00
527 case 3341 : return 0; // SSC* .28000E+01 0.00
528 case -3341 : return 0; // ASSC* .28000E+01 0.00
529 case 1441 : return 0; // UCC* .37500E+01 2.00
530 case -1441 : return 0; // AUCC* .37500E+01 -2.00
531 case 2441 : return 0; // DCC* .37500E+01 1.00
532 case -2441 : return 0; // ADCC* .37500E+01 -1.00
533 case 3441 : return 0; // SCC* .39000E+01 1.00
534 case -3441 : return 0; // ASCC* .39000E+01 -1.00
535 case 4441 : return 0; // CCC* .48000E+01 2.00
536 case -4441 : return 0; // ACCC* .48000E+01 -2.00
537 case 10 : return 22; // Photon
538 case 12 : return 11; // Electron
539 case -12 : return -11; // Positron
540 case 14 : return 13; // Muon-
541 case -14 : return -13; // Muon+
542 case 16 : return 15; // Tau-
543 case -16 : return -15; // Tau+
544 case 11 : return 12; // Neutrino e
545 case -11 : return -12; // Anti Neutrino e
546 case 13 : return 14; // Neutrino Muon
547 case -13 : return -14; // Anti Neutrino Muon
548 case 15 : return 16; // Neutrino Tau
549 case -15 : return -16; // Anti Neutrino Tau
550 case 110 : return 111; // Pion0
551 case 120 : return 211; // Pion+
552 case -120 : return -211; // Pion-
553 case 220 : return 221; // Eta
554 case 130 : return 321; // Kaon+
555 case -130 : return -321; // Kaon-
556 case -20 : return 130; // Kaon Long
557 case 20 : return 310; // Kaon Short
558
559 // baryons
560 case 1120 : return 2212; // Proton
561 case -1120 : return -2212; // Anti Proton
562 case 1220 : return 2112; // Neutron
563 case -1220 : return -2112; // Anti Neutron
564 case 2130 : return 3122; // Lambda
565 case -2130 : return -3122; // Lambda bar
566 case 1130 : return 3222; // Sigma+
567 case -1130 : return -3222; // Sigma bar -
568 case 1230 : return 3212; // Sigma0
569 case -1230 : return -3212; // Sigma bar 0
570 case 2230 : return 3112; // Sigma-
571 case -2230 : return -3112; // Sigma bar +
572 case 1330 : return 3322; // Xi0
573 case -1330 : return -3322; // Xi bar 0
574 case 2330 : return 3312; // Xi-
575 case -2330 : return -3312; // Xi bar +
576 default : return 0; // isajet or pdg number does not exist
577 }
578}
579
581{
582
584 // To avoid to rebuild the particle list from multiple threads, we check if at
585 // the moment of invoked ReadPDGTableImpl, the list was initialized. If it was not, and now it is, another
586 // thread filled it. Doing so again would result in an error.
587 // However, if the particle list was not initialized when ReadPDGTableImpl was invoked and it is still empty,
588 // we fill it.
589 // With the protection described above, ReadPDGTable can be invoked multiple times, for example once
590 // automatically by the system lazily and a second time by the user to complement the information
591 // already in memory with a second input file.
593 return;
594 }
595
596 if (fParticleList == 0) {
599 }
600
602 const char *fn;
603
604 if (!FileName[0]) {
605 default_name = "pdg_table.txt";
607 fn = gEnv->GetValue("Root.DatabasePDG", default_name.Data());
608 } else {
609 fn = FileName;
610 }
611
612 FILE* file = fopen(fn,"r");
613 if (file == 0) {
614 Error("ReadPDGTable","Could not open PDG particle file %s",fn);
615 return;
616 }
617
618 char c[512];
620 Int_t ich, kf, nch, charge;
621 char name[30], class_name[30];
623 Int_t dau[20];
624
626
627 Int_t input;
628 while ( (input=getc(file)) != EOF) {
629 c[0] = input;
630 if (c[0] != '#') {
631 ungetc(c[0],file);
632 // read channel number
633 // coverity [secure_coding : FALSE]
634 if (fscanf(file,"%i",&ich)) {;}
635 // coverity [secure_coding : FALSE]
636 if (fscanf(file,"%s",name )) {;}
637 // coverity [secure_coding : FALSE]
638 if (fscanf(file,"%i",&kf )) {;}
639 // coverity [secure_coding : FALSE]
640 if (fscanf(file,"%i",&anti )) {;}
641
642 if (kf < 0) {
644 // nothing more on this line
645 if (fgets(c,200,file)) {;}
646 } else {
647 // coverity [secure_coding : FALSE]
648 if (fscanf(file,"%i",&class_number)) {;}
649 // coverity [secure_coding : FALSE]
650 if (fscanf(file,"%s",class_name)) {;}
651 // coverity [secure_coding : FALSE]
652 if (fscanf(file,"%i",&charge)) {;}
653 // coverity [secure_coding : FALSE]
654 if (fscanf(file,"%le",&mass)) {;}
655 // coverity [secure_coding : FALSE]
656 if (fscanf(file,"%le",&width)) {;}
657 // coverity [secure_coding : FALSE]
658 if (fscanf(file,"%i",&isospin)) {;}
659 // coverity [secure_coding : FALSE]
660 if (fscanf(file,"%i",&i3)) {;}
661 // coverity [secure_coding : FALSE]
662 if (fscanf(file,"%i",&spin)) {;}
663 // coverity [secure_coding : FALSE]
664 if (fscanf(file,"%i",&flavor)) {;}
665 // coverity [secure_coding : FALSE]
666 if (fscanf(file,"%i",&tracking_code)) {;}
667 // coverity [secure_coding : FALSE]
668 if (fscanf(file,"%i",&nch)) {;}
669 // nothing more on this line
670 if (fgets(c,200,file)) {;}
671 if (width > 1e-10) stable = 0;
672 else stable = 1;
673
674 // create particle
675
677 name,
678 mass,
679 stable,
680 width,
681 charge,
683 kf,
684 anti,
686
687 if (nch) {
688 // read in decay channels
689 ich = 0;
690 Int_t c_input = 0;
691 while ( ((c_input=getc(file)) != EOF) && (ich <nch)) {
692 c[0] = c_input;
693 if (c[0] != '#') {
694 ungetc(c[0],file);
695
696 // coverity [secure_coding : FALSE]
697 if (fscanf(file,"%i",&idecay)) {;}
698 // coverity [secure_coding : FALSE]
699 if (fscanf(file,"%i",&decay_type)) {;}
700 // coverity [secure_coding : FALSE]
701 if (fscanf(file,"%le",&branching_ratio)) {;}
702 // coverity [secure_coding : FALSE]
703 if (fscanf(file,"%i",&ndau)) {;}
704 for (int idau=0; idau<ndau; idau++) {
705 // coverity [secure_coding : FALSE]
706 if (fscanf(file,"%i",&dau[idau])) {;}
707 }
708 // add decay channel
709
710 if (part) part->AddDecayChannel(decay_type,branching_ratio,ndau,dau);
711 ich++;
712 }
713 // skip end of line
714 if (fgets(c,200,file)) {;}
715 }
716 }
717 }
718 } else {
719 // skip end of line
720 if (fgets(c,200,file)) {;}
721 }
722 }
723 // in the end loop over the antiparticles and
724 // define their decay lists
726
727 Int_t code[20];
730
731 while ((p = (TParticlePDG*) it.Next())) {
732
733 // define decay channels for antiparticles
734 if (p->PdgCode() < 0) {
735 ap = GetParticle(-p->PdgCode());
736 if (!ap) continue;
737 nch = ap->NDecayChannels();
738 for (ich=0; ich<nch; ich++) {
739 dc = ap->DecayChannel(ich);
740 if (!dc) continue;
741 ndau = dc->NDaughters();
742 for (int i=0; i<ndau; i++) {
743 // conserve CPT
744
745 code[i] = dc->DaughterPdgCode(i);
746 daughter = GetParticle(code[i]);
747 if (daughter && daughter->AntiParticle()) {
748 // this particle does have an
749 // antiparticle
750 code[i] = -code[i];
751 }
752 }
753 p->AddDecayChannel(dc->MatrixElementCode(),
754 dc->BranchingRatio(),
755 dc->NDaughters(),
756 code);
757 }
758 p->SetAntiParticle(ap);
759 ap->SetAntiParticle(p);
760 }
761 }
762
763 fclose(file);
764 return;
765}
766
767////////////////////////////////////////////////////////////////////////////////
768/// read list of particles from a file
769/// if the particle list does not exist, it is created, otherwise
770/// particles are added to the existing list
771/// See $ROOTSYS/etc/pdg_table.txt to see the file format
772
773void TDatabasePDG::ReadPDGTable(const char *FileName)
774{
775 ReadPDGTableImpl(FileName, nullptr != fParticleList);
776}
777
778////////////////////////////////////////////////////////////////////////////////
779/// browse data base
780
785
786
787////////////////////////////////////////////////////////////////////////////////
788/// write contents of the particle DB into a file
789
791{
792 if (fParticleList == 0) {
793 Error("WritePDGTable","Do not have a valid PDG particle list;"
794 " consider loading it with ReadPDGTable first.");
795 return -1;
796 }
797
798 FILE *file = fopen(filename,"w");
799 if (file == 0) {
800 Error("WritePDGTable","Could not open PDG particle file %s",filename);
801 return -1;
802 }
803
804 fprintf(file,"#--------------------------------------------------------------------\n");
805 fprintf(file,"# i NAME............. KF AP CLASS Q MASS WIDTH 2*I+1 I3 2*S+1 FLVR TrkCod N(dec)\n");
806 fprintf(file,"#--------------------------------------------------------------------\n");
807
809 for(Int_t i=0;i<nparts;++i) {
810 TParticlePDG *p = dynamic_cast<TParticlePDG*>(fParticleList->At(i));
811 if(!p) continue;
812
813 Int_t ich=i+1;
814 Int_t kf=p->PdgCode();
815 fprintf(file,"%5i %-20s %- 6i ", ich, p->GetName(), kf);
816
817 Int_t anti=p->AntiParticle() ? 1:0;
818 if(kf<0) {
819 for(Int_t j=0;j<nparts;++j) {
820 TParticlePDG *dummy = dynamic_cast<TParticlePDG*>(fParticleList->At(j));
821 if(dummy==p->AntiParticle()) {
822 anti=j+1;
823 break;
824 }
825 }
826 fprintf(file,"%i 0\n",anti);
827 continue;
828 }
829
830 fprintf(file,"%i ",anti);
831 fprintf(file,"%i ",100);
832 fprintf(file,"%s ",p->ParticleClass());
833 fprintf(file,"% i ",(Int_t)p->Charge());
834 fprintf(file,"%.5le ",p->Mass());
835 fprintf(file,"%.5le ",p->Width());
836 fprintf(file,"%i ",(Int_t)p->Isospin());
837 fprintf(file,"%i ",(Int_t)p->I3());
838 fprintf(file,"%i ",(Int_t)p->Spin());
839 fprintf(file,"%i ",-1);
840 fprintf(file,"%i ",p->TrackingCode());
841 Int_t nch=p->NDecayChannels();
842 fprintf(file,"%i\n",nch);
843 if(nch==0) {
844 continue;
845 }
846 fprintf(file,"#----------------------------------------------------------------------\n");
847 fprintf(file,"# decay type(PY6) BR Nd daughters(codes, then names)\n");
848 fprintf(file,"#----------------------------------------------------------------------\n");
849 for(Int_t j=0;j<nch;++j) {
850 TDecayChannel *dc=p->DecayChannel(j);
851 if (!dc) continue;
852 fprintf(file,"%9i ",dc->Number()+1);
853 fprintf(file,"%3i ",dc->MatrixElementCode());
854 fprintf(file,"%.5le ",dc->BranchingRatio());
855 Int_t ndau=dc->NDaughters();
856 fprintf(file,"%3i ",ndau);
857 for (int idau=0; idau<ndau; idau++) {
858 fprintf(file,"%- 6i ",dc->DaughterPdgCode(idau));
859 }
860 for (int idau=0; idau<ndau; idau++) {
861 TParticlePDG *dummy=GetParticle(dc->DaughterPdgCode(idau));
862 if(dummy)
863 fprintf(file,"%-10s ",dummy->GetName());
864 else
865 fprintf(file,"%-10s ","???");
866 }
867 fprintf(file,"\n");
868 }
869 }
870 fclose(file);
871 return nparts;
872}
#define b(i)
Definition RSha256.hxx:100
#define c(i)
Definition RSha256.hxx:101
#define e(i)
Definition RSha256.hxx:103
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
TDatabasePDG ** GetInstancePtr()
Static function holding the instance.
R__EXTERN TEnv * gEnv
Definition TEnv.h:170
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t option
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void input
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char filename
Option_t Option_t width
char name[80]
Definition TGX11.cxx:110
#define gROOT
Definition TROOT.h:411
R__EXTERN TSystem * gSystem
Definition TSystem.h:572
#define R__WRITE_LOCKGUARD(mutex)
#define R__READ_LOCKGUARD(mutex)
Using a TBrowser one can browse all ROOT objects.
Definition TBrowser.h:37
virtual Int_t GetEntries() const
void Browse(TBrowser *b) override
Browse this collection (called by TBrowser).
Particle database manager class.
virtual Int_t WritePDGTable(const char *filename)
write contents of the particle DB into a file
static TDatabasePDG * Instance()
static function
virtual Int_t ConvertGeant3ToPdg(Int_t Geant3Number) const
Converts Geant3 particle codes to PDG convention.
virtual Int_t ConvertPdgToGeant3(Int_t pdgNumber) const
Converts pdg code to geant3 id.
void BuildPdgMap() const
Build fPdgMap mapping pdg-code to particle.
virtual Int_t ConvertIsajetToPdg(Int_t isaNumber) const
Converts the ISAJET Particle number into the PDG MC number.
TParticleClassPDG * GetParticleClass(const char *name)
TParticlePDG * GetParticle(Int_t pdgCode) const
Get a pointer to the particle object according to the MC code number.
virtual TParticlePDG * AddParticle(const char *Name, const char *Title, Double_t Mass, Bool_t Stable, Double_t DecayWidth, Double_t Charge, const char *ParticleClass, Int_t PdgCode, Int_t Anti=-1, Int_t TrackingCode=0)
Particle definition normal constructor.
TDatabasePDG()
Create PDG database.
~TDatabasePDG() override
Cleanup the PDG database.
TExMap * fPdgMap
void Print(Option_t *opt="") const override
Print contents of PDG database.
virtual TParticlePDG * AddAntiParticle(const char *Name, Int_t PdgCode)
assuming particle has already been defined
THashList * fParticleList
TObjArray * fListOfClasses
virtual void ReadPDGTableImpl(const char *filename, bool isParticleListInitializedWhenInvoking)
void Browse(TBrowser *b) override
browse data base
virtual void ReadPDGTable(const char *filename="")
read list of particles from a file if the particle list does not exist, it is created,...
Description of the decay channel.
virtual Int_t GetValue(const char *name, Int_t dflt) const
Returns the integer value for a resource.
Definition TEnv.cxx:490
This class stores a (key,value) pair using an external hash.
Definition TExMap.h:33
void Add(ULong64_t hash, Long64_t key, Long64_t value)
Add an (key,value) pair to the table. The key should be unique.
Definition TExMap.cxx:87
Long64_t GetValue(ULong64_t hash, Long64_t key)
Return the value belonging to specified key and hash value.
Definition TExMap.cxx:173
THashList implements a hybrid collection class consisting of a hash table and a list to store TObject...
Definition THashList.h:34
void Delete(Option_t *option="") override
Remove all objects from the list AND delete all heap based objects.
TObject * FindObject(const char *name) const override
Find object using its name.
TObject * Next()
void Add(TObject *obj) override
Definition TList.h:81
TObject * At(Int_t idx) const override
Returns the object at position idx. Returns 0 if idx is out of range.
Definition TList.cxx:354
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
const char * GetName() const override
Returns name of object.
Definition TNamed.h:49
An array of TObjects.
Definition TObjArray.h:31
void Delete(Option_t *option="") override
Remove all objects from the array AND delete all heap based objects.
TObject * Remove(TObject *obj) override
Remove object from array.
void Add(TObject *obj) override
Definition TObjArray.h:68
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:1057
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
@ kInvalidObject
if object ctor succeeded but object should not be used
Definition TObject.h:78
Utility class used internally by TDatabasePDG.
Description of the static properties of a particle.
static const TString & GetEtcDir()
Get the sysconfig directory in the installation. Static utility function.
Definition TROOT.cxx:3107
Basic string class.
Definition TString.h:138
virtual const char * PrependPathName(const char *dir, TString &name)
Concatenate a directory and a file name.
Definition TSystem.cxx:1092
R__EXTERN TVirtualRWMutex * gCoreMutex