Logo ROOT   6.21/01
Reference Guide
TGeoElement.cxx
Go to the documentation of this file.
1 // @(#)root/geom:$Id$
2 // Author: Andrei Gheata 17/06/04
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 /** \class TGeoElement
13 \ingroup Geometry_classes
14 Base class for chemical elements
15 */
16 
17 /** \class TGeoElementRN
18 \ingroup Geometry_classes
19 Class representing a radionuclide
20 */
21 
22 /** \class TGeoElemIter
23 \ingroup Geometry_classes
24 Iterator for decay branches
25 */
26 
27 /** \class TGeoDecayChannel
28 \ingroup Geometry_classes
29 A decay channel for a radionuclide
30 */
31 
32 /** \class TGeoElementTable
33 \ingroup Geometry_classes
34 Table of elements
35 */
36 
37 #include "RConfigure.h"
38 
39 #include "Riostream.h"
40 
41 #include "TSystem.h"
42 #include "TROOT.h"
43 #include "TObjArray.h"
44 #include "TVirtualGeoPainter.h"
45 #include "TGeoManager.h"
46 #include "TGeoElement.h"
47 #include "TMath.h"
48 #include "TGeoPhysicalConstants.h"
50 
51 // statics and globals
52 static const Int_t gMaxElem = 110;
53 static const Int_t gMaxLevel = 8;
54 static const Int_t gMaxDecay = 15;
55 
56 static const char gElName[gMaxElem][3] = {
57  "H ","He","Li","Be","B ","C ","N ","O ","F ","Ne","Na","Mg",
58  "Al","Si","P ","S ","Cl","Ar","K ","Ca","Sc","Ti","V ","Cr",
59  "Mn","Fe","Co","Ni","Cu","Zn","Ga","Ge","As","Se","Br","Kr",
60  "Rb","Sr","Y ","Zr","Nb","Mo","Tc","Ru","Rh","Pd","Ag","Cd",
61  "In","Sn","Sb","Te","I ","Xe","Cs","Ba","La","Ce","Pr","Nd",
62  "Pm","Sm","Eu","Gd","Tb","Dy","Ho","Er","Tm","Yb","Lu","Hf",
63  "Ta","W ","Re","Os","Ir","Pt","Au","Hg","Tl","Pb","Bi","Po",
64  "At","Rn","Fr","Ra","Ac","Th","Pa","U ","Np","Pu","Am","Cm",
65  "Bk","Cf","Es","Fm","Md","No","Lr","Rf","Db","Sg","Bh","Hs",
66  "Mt","Ds" };
67 
68 static const char *gDecayName[gMaxDecay+1] = {
69  "2BetaMinus", "BetaMinus", "NeutronEm", "ProtonEm", "Alpha", "ECF",
70  "ElecCapt", "IsoTrans", "I", "SpontFiss", "2ProtonEm", "2NeutronEm",
71  "2Alpha", "Carbon12", "Carbon14", "Stable" };
72 
73 static const Int_t gDecayDeltaA[gMaxDecay] = {
74  0, 0, -1, -1, -4,
75  -99, 0, 0, -99, -99,
76  -2, -2, -8, -12, -14 };
77 
78 static const Int_t gDecayDeltaZ[gMaxDecay] = {
79  2, 1, 0, -1, -2,
80  -99, -1, 0, -99, -99,
81  -2, 0, -4, -6, -6 };
82 static const char gLevName[gMaxLevel]=" mnpqrs";
83 
85 
86 ////////////////////////////////////////////////////////////////////////////////
87 /// Default constructor
88 
90 {
91  TGeoUnit::setUnitType(TGeoUnit::unitType()); // Ensure nobody changes the units afterwards
93  SetUsed(kFALSE);
94  fZ = 0;
95  fN = 0;
96  fNisotopes = 0;
97  fA = 0.0;
98  fIsotopes = NULL;
99  fAbundances = NULL;
100 }
101 
102 ////////////////////////////////////////////////////////////////////////////////
103 /// Obsolete constructor
104 
105 TGeoElement::TGeoElement(const char *name, const char *title, Int_t z, Double_t a)
106  :TNamed(name, title)
107 {
108  TGeoUnit::setUnitType(TGeoUnit::unitType()); // Ensure nobody changes the units afterwards
110  SetUsed(kFALSE);
111  fZ = z;
112  fN = Int_t(a);
113  fNisotopes = 0;
114  fA = a;
115  fIsotopes = NULL;
116  fAbundances = NULL;
118 }
119 
120 ////////////////////////////////////////////////////////////////////////////////
121 /// Element having isotopes.
122 
123 TGeoElement::TGeoElement(const char *name, const char *title, Int_t nisotopes)
124  :TNamed(name, title)
125 {
126  TGeoUnit::setUnitType(TGeoUnit::unitType()); // Ensure nobody changes the units afterwards
128  SetUsed(kFALSE);
129  fZ = 0;
130  fN = 0;
131  fNisotopes = nisotopes;
132  fA = 0.0;
133  fIsotopes = new TObjArray(nisotopes);
134  fAbundances = new Double_t[nisotopes];
135 }
136 
137 ////////////////////////////////////////////////////////////////////////////////
138 /// Constructor
139 
140 TGeoElement::TGeoElement(const char *name, const char *title, Int_t z, Int_t n, Double_t a)
141  :TNamed(name, title)
142 {
143  TGeoUnit::setUnitType(TGeoUnit::unitType()); // Ensure nobody changes the units afterwards
145  SetUsed(kFALSE);
146  fZ = z;
147  fN = n;
148  fNisotopes = 0;
149  fA = a;
150  fIsotopes = NULL;
151  fAbundances = NULL;
153 }
154 
155 ////////////////////////////////////////////////////////////////////////////////
156 /// destructor
157 
159 {
160  delete fIsotopes;
161  delete [] fAbundances;
162 }
163 
164 
165 ////////////////////////////////////////////////////////////////////////////////
166 /// Calculate properties for an atomic number
167 
169 {
170  // Radiation Length
173 }
174 ////////////////////////////////////////////////////////////////////////////////
175 /// Compute Coulomb correction factor (Phys Rev. D50 3-1 (1994) page 1254)
176 
178 {
179  static constexpr Double_t k1 = 0.0083 , k2 = 0.20206 ,k3 = 0.0020 , k4 = 0.0369;
182  Double_t az2 = (fsc*fZ)*(fsc*fZ);
183  Double_t az4 = az2 * az2;
184 
185  fCoulomb = (k1*az4 + k2 + 1./(1.+az2))*az2 - (k3*az4 + k4)*az4;
186 }
187 ////////////////////////////////////////////////////////////////////////////////
188 /// Compute Tsai's Expression for the Radiation Length (Phys Rev. D50 3-1 (1994) page 1254)
189 
191 {
192  static constexpr Double_t Lrad_light[] = {5.31 , 4.79 , 4.74 , 4.71} ;
193  static constexpr Double_t Lprad_light[] = {6.144 , 5.621 , 5.805 , 5.924} ;
194 
195  fRadTsai = 0.0;
196  if (fZ == 0) return;
197  const Double_t logZ3 = TMath::Log(fZ)/3.;
198 
199  Double_t Lrad, Lprad;
202  Int_t iz = static_cast<Int_t>(fZ+0.5) - 1 ; // The static cast comes from G4lrint
203  static const Double_t log184 = TMath::Log(184.15);
204  static const Double_t log1194 = TMath::Log(1194.);
205  if (iz <= 3) { Lrad = Lrad_light[iz] ; Lprad = Lprad_light[iz] ; }
206  else { Lrad = log184 - logZ3 ; Lprad = log1194 - 2*logZ3;}
207 
208  fRadTsai = 4*alpha_rcl2*fZ*(fZ*(Lrad-fCoulomb) + Lprad);
209 }
210 ////////////////////////////////////////////////////////////////////////////////
211 /// Print this isotope
212 
213 void TGeoElement::Print(Option_t *option) const
214 {
215  printf("Element: %s Z=%d N=%f A=%f [g/mole]\n", GetName(), fZ,Neff(),fA);
216  if (HasIsotopes()) {
217  for (Int_t i=0; i<fNisotopes; i++) {
218  TGeoIsotope *iso = GetIsotope(i);
219  printf("=>Isotope %s, abundance=%f :\n", iso->GetName(), fAbundances[i]);
220  iso->Print(option);
221  }
222  }
223 }
224 
225 ////////////////////////////////////////////////////////////////////////////////
226 /// Returns pointer to the table.
227 
229 {
230  if (!gGeoManager) {
231  ::Error("TGeoElementTable::GetElementTable", "Create a geometry manager first");
232  return NULL;
233  }
234  return gGeoManager->GetElementTable();
235 }
236 
237 ////////////////////////////////////////////////////////////////////////////////
238 /// Add an isotope for this element. All isotopes have to be isotopes of the same element.
239 
240 void TGeoElement::AddIsotope(TGeoIsotope *isotope, Double_t relativeAbundance)
241 {
242  if (!fIsotopes) {
243  Fatal("AddIsotope", "Cannot add isotopes to normal elements. Use constructor with number of isotopes.");
244  return;
245  }
246  Int_t ncurrent = 0;
247  TGeoIsotope *isocrt;
248  for (ncurrent=0; ncurrent<fNisotopes; ncurrent++)
249  if (!fIsotopes->At(ncurrent)) break;
250  if (ncurrent==fNisotopes) {
251  Error("AddIsotope", "All %d isotopes of element %s already defined", fNisotopes, GetName());
252  return;
253  }
254  // Check Z of the new isotope
255  if ((fZ!=0) && (isotope->GetZ()!=fZ)) {
256  Fatal("AddIsotope", "Trying to add isotope %s with different Z to the same element %s",
257  isotope->GetName(), GetName());
258  return;
259  } else {
260  fZ = isotope->GetZ();
261  }
262  fIsotopes->Add(isotope);
263  fAbundances[ncurrent] = relativeAbundance;
264  if (ncurrent==fNisotopes-1) {
265  Double_t weight = 0.0;
266  Double_t aeff = 0.0;
267  Double_t neff = 0.0;
268  for (Int_t i=0; i<fNisotopes; i++) {
269  isocrt = (TGeoIsotope*)fIsotopes->At(i);
270  aeff += fAbundances[i]*isocrt->GetA();
271  neff += fAbundances[i]*isocrt->GetN();
272  weight += fAbundances[i];
273  }
274  aeff /= weight;
275  neff /= weight;
276  fN = (Int_t)neff;
277  fA = aeff;
278  }
280 }
281 
282 ////////////////////////////////////////////////////////////////////////////////
283 /// Returns effective number of nucleons.
284 
286 {
287  if (!fNisotopes) return fN;
288  TGeoIsotope *isocrt;
289  Double_t weight = 0.0;
290  Double_t neff = 0.0;
291  for (Int_t i=0; i<fNisotopes; i++) {
292  isocrt = (TGeoIsotope*)fIsotopes->At(i);
293  neff += fAbundances[i]*isocrt->GetN();
294  weight += fAbundances[i];
295  }
296  neff /= weight;
297  return neff;
298 }
299 
300 ////////////////////////////////////////////////////////////////////////////////
301 /// Return i-th isotope in the element.
302 
304 {
305  if (i>=0 && i<fNisotopes) {
306  return (TGeoIsotope*)fIsotopes->At(i);
307  }
308  return NULL;
309 }
310 
311 ////////////////////////////////////////////////////////////////////////////////
312 /// Return relative abundance of i-th isotope in this element.
313 
315 {
316  if (i>=0 && i<fNisotopes) return fAbundances[i];
317  return 0.0;
318 }
319 
321 
322 ////////////////////////////////////////////////////////////////////////////////
323 /// Dummy I/O constructor
324 
326  :TNamed(),
327  fZ(0),
328  fN(0),
329  fA(0)
330 {
331 }
332 
333 ////////////////////////////////////////////////////////////////////////////////
334 /// Constructor
335 
337  :TNamed(name,""),
338  fZ(z),
339  fN(n),
340  fA(a)
341 {
342  if (z<1) Fatal("ctor", "Not allowed Z=%d (<1) for isotope: %s", z,name);
343  if (n<z) Fatal("ctor", "Not allowed Z=%d < N=%d for isotope: %s", z,n,name);
345 }
346 
347 ////////////////////////////////////////////////////////////////////////////////
348 /// Find existing isotope by name.
349 
351 {
353  if (!elTable) return 0;
354  return elTable->FindIsotope(name);
355 }
356 
357 ////////////////////////////////////////////////////////////////////////////////
358 /// Print this isotope
359 
361 {
362  printf("Isotope: %s Z=%d N=%d A=%f [g/mole]\n", GetName(), fZ,fN,fA);
363 }
364 
366 
367 ////////////////////////////////////////////////////////////////////////////////
368 /// Default constructor
369 
371 {
373  fENDFcode = 0;
374  fIso = 0;
375  fLevel = 0;
376  fDeltaM = 0;
377  fHalfLife = 0;
378  fNatAbun = 0;
379  fTH_F = 0;
380  fTG_F = 0;
381  fTH_S = 0;
382  fTG_S = 0;
383  fStatus = 0;
384  fRatio = 0;
385  fDecays = 0;
386 }
387 
388 ////////////////////////////////////////////////////////////////////////////////
389 /// Constructor.
390 
392  Double_t deltaM, Double_t halfLife, const char* JP,
393  Double_t natAbun, Double_t th_f, Double_t tg_f, Double_t th_s,
394  Double_t tg_s, Int_t status)
395  :TGeoElement("", JP, zz, aa)
396 {
398  fENDFcode = ENDF(aa,zz,iso);
399  fIso = iso;
400  fLevel = level;
401  fDeltaM = deltaM;
402  fHalfLife = halfLife;
403  fTitle = JP;
404  if (!fTitle.Length()) fTitle = "?";
405  fNatAbun = natAbun;
406  fTH_F = th_f;
407  fTG_F = tg_f;
408  fTH_S = th_s;
409  fTG_S = tg_s;
410  fStatus = status;
411  fDecays = 0;
412  fRatio = 0;
413  MakeName(aa,zz,iso);
414  if ((TMath::Abs(fHalfLife)<1.e-30) || fHalfLife<-1) Warning("ctor","Element %s has T1/2=%g [s]", fName.Data(), fHalfLife);
415 }
416 
417 ////////////////////////////////////////////////////////////////////////////////
418 /// Destructor.
419 
421 {
422  if (fDecays) {
423  fDecays->Delete();
424  delete fDecays;
425  }
426  if (fRatio) delete fRatio;
427 }
428 
429 ////////////////////////////////////////////////////////////////////////////////
430 /// Adds a decay mode for this element.
431 
432 void TGeoElementRN::AddDecay(Int_t decay, Int_t diso, Double_t branchingRatio, Double_t qValue)
433 {
434  if (branchingRatio<1E-20) {
435  TString decayName;
436  TGeoDecayChannel::DecayName(decay, decayName);
437  Warning("AddDecay", "Decay %s of %s has BR=0. Not added.", decayName.Data(),fName.Data());
438  return;
439  }
440  TGeoDecayChannel *dc = new TGeoDecayChannel(decay,diso,branchingRatio, qValue);
441  dc->SetParent(this);
442  if (!fDecays) fDecays = new TObjArray(5);
443  fDecays->Add(dc);
444 }
445 
446 ////////////////////////////////////////////////////////////////////////////////
447 /// Adds a decay channel to the list of decays.
448 
450 {
451  dc->SetParent(this);
452  if (!fDecays) fDecays = new TObjArray(5);
453  fDecays->Add(dc);
454 }
455 
456 ////////////////////////////////////////////////////////////////////////////////
457 /// Get number of decay channels of this element.
458 
460 {
461  if (!fDecays) return 0;
462  return fDecays->GetEntriesFast();
463 }
464 
465 ////////////////////////////////////////////////////////////////////////////////
466 /// Get the activity in Bq of a gram of material made from this element.
467 
469 {
470  static const Double_t ln2 = TMath::Log(2.);
471  Double_t sa = (fHalfLife>0 && fA>0)?(ln2*TMath::Na()/fHalfLife/fA):0.;
472  return sa;
473 }
474 
475 ////////////////////////////////////////////////////////////////////////////////
476 /// Check if all decay chain of the element is OK.
477 
479 {
481  TObject *oelem = (TObject*)this;
482  TGeoDecayChannel *dc;
483  TGeoElementRN *elem;
485  TString decayName;
486  if (!table) {
487  Error("CheckDecays", "Element table not present");
488  return kFALSE;
489  }
490  Bool_t resultOK = kTRUE;
491  if (!fDecays) {
492  oelem->SetBit(kElementChecked,kTRUE);
493  return resultOK;
494  }
495  Double_t br = 0.;
496  Int_t decayResult = 0;
497  TIter next(fDecays);
498  while ((dc=(TGeoDecayChannel*)next())) {
499  br += dc->BranchingRatio();
500  decayResult = DecayResult(dc);
501  if (decayResult) {
502  elem = table->GetElementRN(decayResult);
503  if (!elem) {
504  TGeoDecayChannel::DecayName(dc->Decay(),decayName);
505  Error("CheckDecays", "Element after decay %s of %s not found in DB", decayName.Data(),fName.Data());
506  return kFALSE;
507  }
508  dc->SetDaughter(elem);
509  resultOK = elem->CheckDecays();
510  }
511  }
512  if (TMath::Abs(br-100) > 1.E-3) {
513  Warning("CheckDecays", "BR for decays of element %s sum-up = %f", fName.Data(), br);
514  resultOK = kFALSE;
515  }
516  oelem->SetBit(kElementChecked,kTRUE);
517  return resultOK;
518 }
519 
520 ////////////////////////////////////////////////////////////////////////////////
521 /// Returns ENDF code of decay result.
522 
524 {
525  Int_t da, dz, diso;
526  dc->DecayShift(da, dz, diso);
527  if (da == -99 || dz == -99) return 0;
528  return ENDF(Int_t(fA)+da,fZ+dz,fIso+diso);
529 }
530 
531 ////////////////////////////////////////////////////////////////////////////////
532 /// Fills the input array with the set of RN elements resulting from the decay of
533 /// this one. All element in the list will contain the time evolution of their
534 /// proportion by number with respect to this element. The proportion can be
535 /// retrieved via the method TGeoElementRN::Ratio().
536 /// The precision represent the minimum cumulative branching ratio for
537 /// which decay products are still taken into account.
538 
539 void TGeoElementRN::FillPopulation(TObjArray *population, Double_t precision, Double_t factor)
540 {
541  TGeoElementRN *elem;
542  TGeoElemIter next(this, precision);
543  TGeoBatemanSol s(this);
544  s.Normalize(factor);
545  AddRatio(s);
546  if (!population->FindObject(this)) population->Add(this);
547  while ((elem=next())) {
548  TGeoBatemanSol ratio(next.GetBranch());
549  ratio.Normalize(factor);
550  elem->AddRatio(ratio);
551  if (!population->FindObject(elem)) population->Add(elem);
552  }
553 }
554 
555 ////////////////////////////////////////////////////////////////////////////////
556 /// Generate a default name for the element.
557 
559 {
560  fName = "";
561  if (z==0 && a==1) {
562  fName = "neutron";
563  return;
564  }
565  if (z>=1 && z<= gMaxElem) fName += TString::Format("%3d-%s-",z,gElName[z-1]);
566  else fName = "?? -?? -";
567  if (a>=1 && a<=999) fName += TString::Format("%3.3d",a);
568  else fName += "??";
569  if (iso>0 && iso<gMaxLevel) fName += TString::Format("%c", gLevName[iso]);
570  fName.ReplaceAll(" ","");
571 }
572 
573 ////////////////////////////////////////////////////////////////////////////////
574 /// Print info about the element;
575 
576 void TGeoElementRN::Print(Option_t *option) const
577 {
578  printf("\n%-12s ",fName.Data());
579  printf("ENDF=%d; ",fENDFcode);
580  printf("A=%d; ",(Int_t)fA);
581  printf("Z=%d; ",fZ);
582  printf("Iso=%d; ",fIso);
583  printf("Level=%g[MeV]; ",fLevel);
584  printf("Dmass=%g[MeV]; ",fDeltaM);
585  if (fHalfLife>0) printf("Hlife=%g[s]\n",fHalfLife);
586  else printf("Hlife=INF\n");
587  printf("%13s"," ");
588  printf("J/P=%s; ",fTitle.Data());
589  printf("Abund=%g; ",fNatAbun);
590  printf("Htox=%g; ",fTH_F);
591  printf("Itox=%g; ",fTG_F);
592  printf("Stat=%d\n",fStatus);
593  if(!fDecays) return;
594  printf("Decay modes:\n");
595  TIter next(fDecays);
596  TGeoDecayChannel *dc;
597  while ((dc=(TGeoDecayChannel*)next())) dc->Print(option);
598 }
599 
600 ////////////////////////////////////////////////////////////////////////////////
601 /// Create element from line record.
602 
604 {
605  Int_t a,z,iso,status;
606  Double_t level, deltaM, halfLife, natAbun, th_f, tg_f, th_s, tg_s;
607  char name[20],jp[20];
608  sscanf(&line[0], "%s%d%d%d%lg%lg%lg%s%lg%lg%lg%lg%lg%d%d", name,&a,&z,&iso,&level,&deltaM,
609  &halfLife,jp,&natAbun,&th_f,&tg_f,&th_s,&tg_s,&status,&ndecays);
610  TGeoElementRN *elem = new TGeoElementRN(a,z,iso,level,deltaM,halfLife,
611  jp,natAbun,th_f,tg_f,th_s,tg_s,status);
612  return elem;
613 }
614 
615 ////////////////////////////////////////////////////////////////////////////////
616 /// Save primitive for RN elements.
617 
618 void TGeoElementRN::SavePrimitive(std::ostream &out, Option_t *option)
619 {
620  if (!strcmp(option,"h")) {
621  // print a header if requested
622  out << "#====================================================================================================================================" << std::endl;
623  out << "# Name A Z ISO LEV[MeV] DM[MeV] T1/2[s] J/P ABUND[%] HTOX ITOX HTOX ITOX STAT NDCY" << std::endl;
624  out << "#====================================================================================================================================" << std::endl;
625  }
626  out << std::setw(11) << fName.Data();
627  out << std::setw(5) << (Int_t)fA;
628  out << std::setw(5) << fZ;
629  out << std::setw(5) << fIso;
630  out << std::setw(10) << std::setiosflags(std::ios::fixed) << std::setprecision(5) << fLevel;
631  out << std::setw(10) << std::setiosflags(std::ios::fixed) << std::setprecision(5) << fDeltaM;
632  out << std::setw(10) << std::setiosflags(std::ios::scientific) << std::setprecision(3) << fHalfLife;
633  out << std::setw(13) << fTitle.Data();
634  out << std::setw(10) << std::setiosflags(std::ios::fixed) << std::setprecision(5) << fNatAbun;
635  out << std::setw(10) << std::setiosflags(std::ios::fixed) << std::setprecision(5) << fTH_F;
636  out << std::setw(10) << std::setiosflags(std::ios::fixed) << std::setprecision(5) << fTG_F;
637  out << std::setw(10) << std::setiosflags(std::ios::fixed) << std::setprecision(5) << fTH_S;
638  out << std::setw(10) << std::setiosflags(std::ios::fixed) << std::setprecision(5) << fTG_S;
639  out << std::setw(5) << fStatus;
640  Int_t ndecays = 0;
641  if (fDecays) ndecays = fDecays->GetEntries();
642  out << std::setw(5) << ndecays;
643  out << std::endl;
644  if (fDecays) {
645  TIter next(fDecays);
646  TGeoDecayChannel *dc;
647  while ((dc=(TGeoDecayChannel*)next())) dc->SavePrimitive(out);
648  }
649 }
650 
651 ////////////////////////////////////////////////////////////////////////////////
652 /// Adds a proportion ratio to the existing one.
653 
655 {
656  if (!fRatio) fRatio = new TGeoBatemanSol(ratio);
657  else *fRatio += ratio;
658 }
659 
660 ////////////////////////////////////////////////////////////////////////////////
661 /// Clears the existing ratio.
662 
664 {
665  if (fRatio) {
666  delete fRatio;
667  fRatio = 0;
668  }
669 }
670 
672 
673 ////////////////////////////////////////////////////////////////////////////////
674 /// Assignment.
675 ///assignment operator
676 
678 {
679  if(this!=&dc) {
680  TObject::operator=(dc);
681  fDecay = dc.fDecay;
682  fDiso = dc.fDiso;
684  fQvalue = dc.fQvalue;
685  fParent = dc.fParent;
686  fDaughter = dc.fDaughter;
687  }
688  return *this;
689 }
690 
691 ////////////////////////////////////////////////////////////////////////////////
692 /// Returns name of decay.
693 
694 const char *TGeoDecayChannel::GetName() const
695 {
696  static TString name = "";
697  name = "";
698  if (!fDecay) return gDecayName[gMaxDecay];
699  for (Int_t i=0; i<gMaxDecay; i++) {
700  if (1<<i & fDecay) {
701  if (name.Length()) name += "+";
702  name += gDecayName[i];
703  }
704  }
705  return name.Data();
706 }
707 
708 ////////////////////////////////////////////////////////////////////////////////
709 /// Returns decay name.
710 
712 {
713  if (!decay) {
715  return;
716  }
717  name = "";
718  for (Int_t i=0; i<gMaxDecay; i++) {
719  if (1<<i & decay) {
720  if (name.Length()) name += "+";
721  name += gDecayName[i];
722  }
723  }
724 }
725 
726 ////////////////////////////////////////////////////////////////////////////////
727 /// Get index of this channel in the list of decays of the parent nuclide.
728 
730 {
731  return fParent->Decays()->IndexOf(this);
732 }
733 
734 ////////////////////////////////////////////////////////////////////////////////
735 /// Prints decay info.
736 
738 {
739  TString name;
741  printf("%-20s Diso: %3d BR: %9.3f%% Qval: %g\n", name.Data(),fDiso,fBranchingRatio,fQvalue);
742 }
743 
744 ////////////////////////////////////////////////////////////////////////////////
745 /// Create element from line record.
746 
748 {
749  char name[80];
750  Int_t decay,diso;
751  Double_t branchingRatio, qValue;
752  sscanf(&line[0], "%s%d%d%lg%lg", name,&decay,&diso,&branchingRatio,&qValue);
753  TGeoDecayChannel *dc = new TGeoDecayChannel(decay,diso,branchingRatio,qValue);
754  return dc;
755 }
756 
757 ////////////////////////////////////////////////////////////////////////////////
758 /// Save primitive for decays.
759 
760 void TGeoDecayChannel::SavePrimitive(std::ostream &out, Option_t *)
761 {
762  TString decayName;
763  DecayName(fDecay, decayName);
764  out << std::setw(50) << decayName.Data();
765  out << std::setw(10) << fDecay;
766  out << std::setw(10) << fDiso;
767  out << std::setw(12) << std::setiosflags(std::ios::fixed) << std::setprecision(6) << fBranchingRatio;
768  out << std::setw(12) << std::setiosflags(std::ios::fixed) << std::setprecision(6) << fQvalue;
769  out << std::endl;
770 }
771 
772 ////////////////////////////////////////////////////////////////////////////////
773 /// Returns variation in A, Z and Iso after decay.
774 
776 {
777  dA=dZ=0;
778  dI=fDiso;
779  for(Int_t i=0; i<gMaxDecay; ++i) {
780  if(1<<i & fDecay) {
781  if(gDecayDeltaA[i] == -99 || gDecayDeltaZ[i] == -99 ) {
782  dA=dZ=-99;
783  return;
784  }
785  dA += gDecayDeltaA[i];
786  dZ += gDecayDeltaZ[i];
787  }
788  }
789 }
790 
792 
793 ////////////////////////////////////////////////////////////////////////////////
794 /// Default constructor.
795 
797  : fTop(top), fElem(top), fBranch(0), fLevel(0), fLimitRatio(limit), fRatio(1.)
798 {
799  fBranch = new TObjArray(10);
800 }
801 
802 ////////////////////////////////////////////////////////////////////////////////
803 /// Copy ctor.
804 
806  :fTop(iter.fTop),
807  fElem(iter.fElem),
808  fBranch(0),
809  fLevel(iter.fLevel),
810  fLimitRatio(iter.fLimitRatio),
811  fRatio(iter.fRatio)
812 {
813  if (iter.fBranch) {
814  fBranch = new TObjArray(10);
815  for (Int_t i=0; i<fLevel; i++) fBranch->Add(iter.fBranch->At(i));
816  }
817 }
818 
819 ////////////////////////////////////////////////////////////////////////////////
820 /// Destructor.
821 
823 {
824  if (fBranch) delete fBranch;
825 }
826 
827 ////////////////////////////////////////////////////////////////////////////////
828 /// Assignment.
829 
831 {
832  if (&iter == this) return *this;
833  fTop = iter.fTop;
834  fElem = iter.fElem;
835  fLevel = iter.fLevel;
836  if (iter.fBranch) {
837  fBranch = new TObjArray(10);
838  for (Int_t i=0; i<fLevel; i++) fBranch->Add(iter.fBranch->At(i));
839  }
840  fLimitRatio = iter.fLimitRatio;
841  fRatio = iter.fRatio;
842  return *this;
843 }
844 
845 ////////////////////////////////////////////////////////////////////////////////
846 /// () operator.
847 
849 {
850  return Next();
851 }
852 
853 ////////////////////////////////////////////////////////////////////////////////
854 /// Go upwards from the current location until the next branching, then down.
855 
857 {
858  TGeoDecayChannel *dc;
859  Int_t ind, nd;
860  while (fLevel) {
861  // Current decay channel
862  dc = (TGeoDecayChannel*)fBranch->At(fLevel-1);
863  ind = dc->GetIndex();
864  nd = dc->Parent()->GetNdecays();
865  fRatio /= 0.01*dc->BranchingRatio();
866  fElem = dc->Parent();
867  fBranch->RemoveAt(--fLevel);
868  ind++;
869  while (ind<nd) {
870  if (Down(ind++)) return (TGeoElementRN*)fElem;
871  }
872  }
873  fElem = NULL;
874  return NULL;
875 }
876 
877 ////////////////////////////////////////////////////////////////////////////////
878 /// Go downwards from current level via ibranch as low in the tree as possible.
879 /// Return value flags if the operation was successful.
880 
882 {
883  TGeoDecayChannel *dc = (TGeoDecayChannel*)fElem->Decays()->At(ibranch);
884  if (!dc->Daughter()) return NULL;
885  Double_t br = 0.01*fRatio*dc->BranchingRatio();
886  if (br < fLimitRatio) return NULL;
887  fLevel++;
888  fRatio = br;
889  fBranch->Add(dc);
890  fElem = dc->Daughter();
891  return (TGeoElementRN*)fElem;
892 }
893 
894 ////////////////////////////////////////////////////////////////////////////////
895 /// Return next element.
896 
898 {
899  if (!fElem) return NULL;
900  // Check if this is the first iteration.
901  Int_t nd = fElem->GetNdecays();
902  for (Int_t i=0; i<nd; i++) if (Down(i)) return (TGeoElementRN*)fElem;
903  return Up();
904 }
905 
906 ////////////////////////////////////////////////////////////////////////////////
907 /// Print info about the current decay branch.
908 
909 void TGeoElemIter::Print(Option_t * /*option*/) const
910 {
911  TGeoElementRN *elem;
912  TGeoDecayChannel *dc;
913  TString indent = "";
914  printf("=== Chain with %g %%\n", 100*fRatio);
915  for (Int_t i=0; i<fLevel; i++) {
916  dc = (TGeoDecayChannel*)fBranch->At(i);
917  elem = dc->Parent();
918  printf("%s%s (%g%% %s) T1/2=%g\n", indent.Data(), elem->GetName(),dc->BranchingRatio(),dc->GetName(),elem->HalfLife());
919  indent += " ";
920  if (i==fLevel-1) {
921  elem = dc->Daughter();
922  printf("%s%s\n", indent.Data(), elem->GetName());
923  }
924  }
925 }
926 
928 
929 ////////////////////////////////////////////////////////////////////////////////
930 /// default constructor
931 
933 {
934  fNelements = 0;
935  fNelementsRN = 0;
936  fNisotopes = 0;
937  fList = 0;
938  fListRN = 0;
939  fIsotopes = 0;
940 }
941 
942 ////////////////////////////////////////////////////////////////////////////////
943 /// constructor
944 
946 {
947  fNelements = 0;
948  fNelementsRN = 0;
949  fNisotopes = 0;
950  fList = new TObjArray(128);
951  fListRN = 0;
952  fIsotopes = 0;
954 // BuildElementsRN();
955 }
956 
957 ////////////////////////////////////////////////////////////////////////////////
958 ///copy constructor
959 
961  TObject(get),
962  fNelements(get.fNelements),
963  fNelementsRN(get.fNelementsRN),
964  fNisotopes(get.fNisotopes),
965  fList(get.fList),
966  fListRN(get.fListRN),
967  fIsotopes(0)
968 {
969 }
970 
971 ////////////////////////////////////////////////////////////////////////////////
972 ///assignment operator
973 
975 {
976  if(this!=&get) {
977  TObject::operator=(get);
978  fNelements=get.fNelements;
979  fNelementsRN=get.fNelementsRN;
980  fNisotopes=get.fNisotopes;
981  fList=get.fList;
982  fListRN=get.fListRN;
983  fIsotopes = 0;
984  }
985  return *this;
986 }
987 
988 ////////////////////////////////////////////////////////////////////////////////
989 /// destructor
990 
992 {
993  if (fList) {
994  fList->Delete();
995  delete fList;
996  }
997  if (fListRN) {
998  fListRN->Delete();
999  delete fListRN;
1000  }
1001  if (fIsotopes) {
1002  fIsotopes->Delete();
1003  delete fIsotopes;
1004  }
1005 }
1006 
1007 ////////////////////////////////////////////////////////////////////////////////
1008 /// Add an element to the table. Obsolete.
1009 
1010 void TGeoElementTable::AddElement(const char *name, const char *title, Int_t z, Double_t a)
1011 {
1012  if (!fList) fList = new TObjArray(128);
1013  fList->AddAtAndExpand(new TGeoElement(name,title,z,a), fNelements++);
1014 }
1015 
1016 ////////////////////////////////////////////////////////////////////////////////
1017 /// Add an element to the table.
1018 
1019 void TGeoElementTable::AddElement(const char *name, const char *title, Int_t z, Int_t n, Double_t a)
1020 {
1021  if (!fList) fList = new TObjArray(128);
1022  fList->AddAtAndExpand(new TGeoElement(name,title,z,n,a), fNelements++);
1023 }
1024 
1025 ////////////////////////////////////////////////////////////////////////////////
1026 /// Add a custom element to the table.
1027 
1029 {
1030  if (!fList) fList = new TObjArray(128);
1031  TGeoElement *orig = FindElement(elem->GetName());
1032  if (orig) {
1033  Error("AddElement", "Found element with same name: %s (%s). Cannot add to table.",
1034  orig->GetName(), orig->GetTitle());
1035  return;
1036  }
1037  fList->AddAtAndExpand(elem, fNelements++);
1038 }
1039 
1040 ////////////////////////////////////////////////////////////////////////////////
1041 /// Add a radionuclide to the table and map it.
1042 
1044 {
1045  if (!fListRN) fListRN = new TObjArray(3600);
1046  if (HasRNElements() && GetElementRN(elem->ENDFCode())) return;
1047 // elem->Print();
1048  fListRN->Add(elem);
1049  fNelementsRN++;
1050  fElementsRN.insert(ElementRNMap_t::value_type(elem->ENDFCode(), elem));
1051 }
1052 
1053 ////////////////////////////////////////////////////////////////////////////////
1054 /// Add isotope to the table.
1055 
1057 {
1058  if (FindIsotope(isotope->GetName())) {
1059  Error("AddIsotope", "Isotope with the same name: %s already in table. Not adding.",isotope->GetName());
1060  return;
1061  }
1062  if (!fIsotopes) fIsotopes = new TObjArray();
1063  fIsotopes->Add(isotope);
1064 }
1065 
1066 ////////////////////////////////////////////////////////////////////////////////
1067 /// Creates the default element table
1068 
1070 {
1071  if (HasDefaultElements()) return;
1072  AddElement("VACUUM","VACUUM" ,0, 0, 0.0);
1073  AddElement("H" ,"HYDROGEN" ,1, 1, 1.00794);
1074  AddElement("HE" ,"HELIUM" ,2, 4, 4.002602);
1075  AddElement("LI" ,"LITHIUM" ,3, 7, 6.941);
1076  AddElement("BE" ,"BERYLLIUM" ,4, 9, 9.01218);
1077  AddElement("B" ,"BORON" ,5, 11, 10.811);
1078  AddElement("C" ,"CARBON" ,6, 12, 12.0107);
1079  AddElement("N" ,"NITROGEN" ,7, 14, 14.00674);
1080  AddElement("O" ,"OXYGEN" ,8, 16, 15.9994);
1081  AddElement("F" ,"FLUORINE" ,9, 19, 18.9984032);
1082  AddElement("NE" ,"NEON" ,10, 20, 20.1797);
1083  AddElement("NA" ,"SODIUM" ,11, 23, 22.989770);
1084  AddElement("MG" ,"MAGNESIUM" ,12, 24, 24.3050);
1085  AddElement("AL" ,"ALUMINIUM" ,13, 27, 26.981538);
1086  AddElement("SI" ,"SILICON" ,14, 28, 28.0855);
1087  AddElement("P" ,"PHOSPHORUS" ,15, 31, 30.973761);
1088  AddElement("S" ,"SULFUR" ,16, 32, 32.066);
1089  AddElement("CL" ,"CHLORINE" ,17, 35, 35.4527);
1090  AddElement("AR" ,"ARGON" ,18, 40, 39.948);
1091  AddElement("K" ,"POTASSIUM" ,19, 39, 39.0983);
1092  AddElement("CA" ,"CALCIUM" ,20, 40, 40.078);
1093  AddElement("SC" ,"SCANDIUM" ,21, 45, 44.955910);
1094  AddElement("TI" ,"TITANIUM" ,22, 48, 47.867);
1095  AddElement("V" ,"VANADIUM" ,23, 51, 50.9415);
1096  AddElement("CR" ,"CHROMIUM" ,24, 52, 51.9961);
1097  AddElement("MN" ,"MANGANESE" ,25, 55, 54.938049);
1098  AddElement("FE" ,"IRON" ,26, 56, 55.845);
1099  AddElement("CO" ,"COBALT" ,27, 59, 58.933200);
1100  AddElement("NI" ,"NICKEL" ,28, 59, 58.6934);
1101  AddElement("CU" ,"COPPER" ,29, 64, 63.546);
1102  AddElement("ZN" ,"ZINC" ,30, 65, 65.39);
1103  AddElement("GA" ,"GALLIUM" ,31, 70, 69.723);
1104  AddElement("GE" ,"GERMANIUM" ,32, 73, 72.61);
1105  AddElement("AS" ,"ARSENIC" ,33, 75, 74.92160);
1106  AddElement("SE" ,"SELENIUM" ,34, 79, 78.96);
1107  AddElement("BR" ,"BROMINE" ,35, 80, 79.904);
1108  AddElement("KR" ,"KRYPTON" ,36, 84, 83.80);
1109  AddElement("RB" ,"RUBIDIUM" ,37, 85, 85.4678);
1110  AddElement("SR" ,"STRONTIUM" ,38, 88, 87.62);
1111  AddElement("Y" ,"YTTRIUM" ,39, 89, 88.90585);
1112  AddElement("ZR" ,"ZIRCONIUM" ,40, 91, 91.224);
1113  AddElement("NB" ,"NIOBIUM" ,41, 93, 92.90638);
1114  AddElement("MO" ,"MOLYBDENUM" ,42, 96, 95.94);
1115  AddElement("TC" ,"TECHNETIUM" ,43, 98, 98.0);
1116  AddElement("RU" ,"RUTHENIUM" ,44, 101, 101.07);
1117  AddElement("RH" ,"RHODIUM" ,45, 103, 102.90550);
1118  AddElement("PD" ,"PALLADIUM" ,46, 106, 106.42);
1119  AddElement("AG" ,"SILVER" ,47, 108, 107.8682);
1120  AddElement("CD" ,"CADMIUM" ,48, 112, 112.411);
1121  AddElement("IN" ,"INDIUM" ,49, 115, 114.818);
1122  AddElement("SN" ,"TIN" ,50, 119, 118.710);
1123  AddElement("SB" ,"ANTIMONY" ,51, 122, 121.760);
1124  AddElement("TE" ,"TELLURIUM" ,52, 128, 127.60);
1125  AddElement("I" ,"IODINE" ,53, 127, 126.90447);
1126  AddElement("XE" ,"XENON" ,54, 131, 131.29);
1127  AddElement("CS" ,"CESIUM" ,55, 133, 132.90545);
1128  AddElement("BA" ,"BARIUM" ,56, 137, 137.327);
1129  AddElement("LA" ,"LANTHANUM" ,57, 139, 138.9055);
1130  AddElement("CE" ,"CERIUM" ,58, 140, 140.116);
1131  AddElement("PR" ,"PRASEODYMIUM" ,59, 141, 140.90765);
1132  AddElement("ND" ,"NEODYMIUM" ,60, 144, 144.24);
1133  AddElement("PM" ,"PROMETHIUM" ,61, 145, 145.0);
1134  AddElement("SM" ,"SAMARIUM" ,62, 150, 150.36);
1135  AddElement("EU" ,"EUROPIUM" ,63, 152, 151.964);
1136  AddElement("GD" ,"GADOLINIUM" ,64, 157, 157.25);
1137  AddElement("TB" ,"TERBIUM" ,65, 159, 158.92534);
1138  AddElement("DY" ,"DYSPROSIUM" ,66, 162, 162.50);
1139  AddElement("HO" ,"HOLMIUM" ,67, 165, 164.93032);
1140  AddElement("ER" ,"ERBIUM" ,68, 167, 167.26);
1141  AddElement("TM" ,"THULIUM" ,69, 169, 168.93421);
1142  AddElement("YB" ,"YTTERBIUM" ,70, 173, 173.04);
1143  AddElement("LU" ,"LUTETIUM" ,71, 175, 174.967);
1144  AddElement("HF" ,"HAFNIUM" ,72, 178, 178.49);
1145  AddElement("TA" ,"TANTALUM" ,73, 181, 180.9479);
1146  AddElement("W" ,"TUNGSTEN" ,74, 184, 183.84);
1147  AddElement("RE" ,"RHENIUM" ,75, 186, 186.207);
1148  AddElement("OS" ,"OSMIUM" ,76, 190, 190.23);
1149  AddElement("IR" ,"IRIDIUM" ,77, 192, 192.217);
1150  AddElement("PT" ,"PLATINUM" ,78, 195, 195.078);
1151  AddElement("AU" ,"GOLD" ,79, 197, 196.96655);
1152  AddElement("HG" ,"MERCURY" ,80, 200, 200.59);
1153  AddElement("TL" ,"THALLIUM" ,81, 204, 204.3833);
1154  AddElement("PB" ,"LEAD" ,82, 207, 207.2);
1155  AddElement("BI" ,"BISMUTH" ,83, 209, 208.98038);
1156  AddElement("PO" ,"POLONIUM" ,84, 209, 209.0);
1157  AddElement("AT" ,"ASTATINE" ,85, 210, 210.0);
1158  AddElement("RN" ,"RADON" ,86, 222, 222.0);
1159  AddElement("FR" ,"FRANCIUM" ,87, 223, 223.0);
1160  AddElement("RA" ,"RADIUM" ,88, 226, 226.0);
1161  AddElement("AC" ,"ACTINIUM" ,89, 227, 227.0);
1162  AddElement("TH" ,"THORIUM" ,90, 232, 232.0381);
1163  AddElement("PA" ,"PROTACTINIUM" ,91, 231, 231.03588);
1164  AddElement("U" ,"URANIUM" ,92, 238, 238.0289);
1165  AddElement("NP" ,"NEPTUNIUM" ,93, 237, 237.0);
1166  AddElement("PU" ,"PLUTONIUM" ,94, 244, 244.0);
1167  AddElement("AM" ,"AMERICIUM" ,95, 243, 243.0);
1168  AddElement("CM" ,"CURIUM" ,96, 247, 247.0);
1169  AddElement("BK" ,"BERKELIUM" ,97, 247, 247.0);
1170  AddElement("CF" ,"CALIFORNIUM",98, 251, 251.0);
1171  AddElement("ES" ,"EINSTEINIUM",99, 252, 252.0);
1172  AddElement("FM" ,"FERMIUM" ,100, 257, 257.0);
1173  AddElement("MD" ,"MENDELEVIUM",101, 258, 258.0);
1174  AddElement("NO" ,"NOBELIUM" ,102, 259, 259.0);
1175  AddElement("LR" ,"LAWRENCIUM" ,103, 262, 262.0);
1176  AddElement("RF" ,"RUTHERFORDIUM",104, 261, 261.0);
1177  AddElement("DB" ,"DUBNIUM" ,105, 262, 262.0);
1178  AddElement("SG" ,"SEABORGIUM" ,106, 263, 263.0);
1179  AddElement("BH" ,"BOHRIUM" ,107, 262, 262.0);
1180  AddElement("HS" ,"HASSIUM" ,108, 265, 265.0);
1181  AddElement("MT" ,"MEITNERIUM" ,109, 266, 266.0);
1182  AddElement("UUN" ,"UNUNNILIUM" ,110, 269, 269.0);
1183  AddElement("UUU" ,"UNUNUNIUM" ,111, 272, 272.0);
1184  AddElement("UUB" ,"UNUNBIUM" ,112, 277, 277.0);
1185 
1187 }
1188 
1189 ////////////////////////////////////////////////////////////////////////////////
1190 /// Creates the list of radionuclides.
1191 
1193 {
1194  if (HasRNElements()) return;
1195  TGeoElementRN *elem;
1196  TString rnf = "RadioNuclides.txt";
1198  FILE *fp = fopen(rnf, "r");
1199  if (!fp) {
1200  Error("ImportElementsRN","File RadioNuclides.txt not found");
1201  return;
1202  }
1203  char line[150];
1204  Int_t ndecays = 0;
1205  Int_t i;
1206  while (fgets(&line[0],140,fp)) {
1207  if (line[0]=='#') continue;
1208  elem = TGeoElementRN::ReadElementRN(line, ndecays);
1209  for (i=0; i<ndecays; i++) {
1210  if (!fgets(&line[0],140,fp)) {
1211  Error("ImportElementsRN", "Error parsing RadioNuclides.txt file");
1212  fclose(fp);
1213  return;
1214  }
1216  elem->AddDecay(dc);
1217  }
1218  AddElementRN(elem);
1219 // elem->Print();
1220  }
1222  CheckTable();
1223  fclose(fp);
1224 }
1225 
1226 ////////////////////////////////////////////////////////////////////////////////
1227 /// Checks status of element table.
1228 
1230 {
1231  if (!HasRNElements()) return HasDefaultElements();
1232  TGeoElementRN *elem;
1233  Bool_t result = kTRUE;
1234  TIter next(fListRN);
1235  while ((elem=(TGeoElementRN*)next())) {
1236  if (!elem->CheckDecays()) result = kFALSE;
1237  }
1238  return result;
1239 }
1240 
1241 ////////////////////////////////////////////////////////////////////////////////
1242 /// Export radionuclides in a file.
1243 
1244 void TGeoElementTable::ExportElementsRN(const char *filename)
1245 {
1246  if (!HasRNElements()) return;
1247  TString sname = filename;
1248  if (!sname.Length()) sname = "RadioNuclides.txt";
1249  std::ofstream out;
1250  out.open(sname.Data(), std::ios::out);
1251  if (!out.good()) {
1252  Error("ExportElementsRN", "Cannot open file %s", sname.Data());
1253  return;
1254  }
1255 
1256  TGeoElementRN *elem;
1257  TIter next(fListRN);
1258  Int_t i=0;
1259  while ((elem=(TGeoElementRN*)next())) {
1260  if ((i%48)==0) elem->SavePrimitive(out,"h");
1261  else elem->SavePrimitive(out);
1262  i++;
1263  }
1264  out.close();
1265 }
1266 
1267 ////////////////////////////////////////////////////////////////////////////////
1268 /// Search an element by symbol or full name
1269 /// Exact matching
1270 
1272 {
1273  TGeoElement *elem;
1274  elem = (TGeoElement*)fList->FindObject(name);
1275  if (elem) return elem;
1276  // Search case insensitive by element name
1277  TString s(name);
1278  s.ToUpper();
1279  elem = (TGeoElement*)fList->FindObject(s.Data());
1280  if (elem) return elem;
1281  // Search by full name
1282  TIter next(fList);
1283  while ((elem=(TGeoElement*)next())) {
1284  if (s == elem->GetTitle()) return elem;
1285  }
1286  return 0;
1287 }
1288 
1289 ////////////////////////////////////////////////////////////////////////////////
1290 /// Find existing isotope by name. Not optimized for a big number of isotopes.
1291 
1293 {
1294  if (!fIsotopes) return NULL;
1295  return (TGeoIsotope*)fIsotopes->FindObject(name);
1296 }
1297 
1298 ////////////////////////////////////////////////////////////////////////////////
1299 /// Retrieve a radionuclide by ENDF code.
1300 
1302 {
1303  if (!HasRNElements()) {
1304  TGeoElementTable *table = (TGeoElementTable*)this;
1305  table->ImportElementsRN();
1306  if (!fListRN) return 0;
1307  }
1308  ElementRNMap_t::const_iterator it = fElementsRN.find(ENDFcode);
1309  if (it != fElementsRN.end()) return it->second;
1310  return 0;
1311 }
1312 
1313 ////////////////////////////////////////////////////////////////////////////////
1314 /// Retrieve a radionuclide by a, z, and isomeric state.
1315 
1317 {
1318  return GetElementRN(TGeoElementRN::ENDF(a,z,iso));
1319 }
1320 
1321 ////////////////////////////////////////////////////////////////////////////////
1322 /// Print table of elements. The accepted options are:
1323 /// "" - prints everything by default
1324 /// "D" - prints default elements only
1325 /// "I" - prints isotopes
1326 /// "R" - prints radio-nuclides only if imported
1327 /// "U" - prints user-defined elements only
1328 
1330 {
1331  TString opt(option);
1332  opt.ToUpper();
1333  Int_t induser = HasDefaultElements() ? 113 : 0;
1334  // Default elements
1335  if (opt=="" || opt=="D") {
1336  if (induser) printf("================\nDefault elements\n================\n");
1337  for (Int_t iel=0; iel<induser; ++iel) fList->At(iel)->Print();
1338  }
1339  // Isotopes
1340  if (opt=="" || opt=="I") {
1341  if (fIsotopes) {
1342  printf("================\nIsotopes\n================\n");
1343  fIsotopes->Print();
1344  }
1345  }
1346  // Radio-nuclides
1347  if (opt=="" || opt=="R") {
1348  if (HasRNElements()) {
1349  printf("================\nRadio-nuclides\n================\n");
1350  fListRN->Print();
1351  }
1352  }
1353  // User-defined elements
1354  if (opt=="" || opt=="U") {
1355  if (fNelements>induser) printf("================\nUser elements\n================\n");
1356  for (Int_t iel=induser; iel<fNelements; ++iel) fList->At(iel)->Print();
1357  }
1358 }
1359 
1361 
1362 ////////////////////////////////////////////////////////////////////////////////
1363 /// Default ctor.
1364 
1366  :TObject(), TAttLine(), TAttFill(), TAttMarker(),
1367  fElem(elem),
1368  fElemTop(elem),
1369  fCsize(10),
1370  fNcoeff(0),
1371  fFactor(1.),
1372  fTmin(0.),
1373  fTmax(0.),
1374  fCoeff(NULL)
1375 {
1376  fCoeff = new BtCoef_t[fCsize];
1377  fNcoeff = 1;
1378  fCoeff[0].cn = 1.;
1379  Double_t t12 = elem->HalfLife();
1380  if (t12 == 0.) t12 = 1.e-30;
1381  if (elem->Stable()) fCoeff[0].lambda = 0.;
1382  else fCoeff[0].lambda = TMath::Log(2.)/t12;
1383 }
1384 
1385 ////////////////////////////////////////////////////////////////////////////////
1386 /// Default ctor.
1387 
1389  :TObject(), TAttLine(), TAttFill(), TAttMarker(),
1390  fElem(NULL),
1391  fElemTop(NULL),
1392  fCsize(0),
1393  fNcoeff(0),
1394  fFactor(1.),
1395  fTmin(0.),
1396  fTmax(0.),
1397  fCoeff(NULL)
1398 {
1399  TGeoDecayChannel *dc = (TGeoDecayChannel*)chain->At(0);
1400  if (dc) fElemTop = dc->Parent();
1401  dc = (TGeoDecayChannel*)chain->At(chain->GetEntriesFast()-1);
1402  if (dc) {
1403  fElem = dc->Daughter();
1404  fCsize = chain->GetEntriesFast()+1;
1405  fCoeff = new BtCoef_t[fCsize];
1406  FindSolution(chain);
1407  }
1408 }
1409 
1410 ////////////////////////////////////////////////////////////////////////////////
1411 /// Copy constructor.
1412 
1414  :TObject(other), TAttLine(other), TAttFill(other), TAttMarker(other),
1415  fElem(other.fElem),
1416  fElemTop(other.fElemTop),
1417  fCsize(other.fCsize),
1418  fNcoeff(other.fNcoeff),
1419  fFactor(other.fFactor),
1420  fTmin(other.fTmin),
1421  fTmax(other.fTmax),
1422  fCoeff(NULL)
1423 {
1424  if (fCsize) {
1425  fCoeff = new BtCoef_t[fCsize];
1426  for (Int_t i=0; i<fNcoeff; i++) {
1427  fCoeff[i].cn = other.fCoeff[i].cn;
1428  fCoeff[i].lambda = other.fCoeff[i].lambda;
1429  }
1430  }
1431 }
1432 
1433 ////////////////////////////////////////////////////////////////////////////////
1434 /// Destructor.
1435 
1437 {
1438  if (fCoeff) delete [] fCoeff;
1439 }
1440 
1441 ////////////////////////////////////////////////////////////////////////////////
1442 /// Assignment.
1443 
1445 {
1446  if (this == &other) return *this;
1447  TObject::operator=(other);
1448  TAttLine::operator=(other);
1449  TAttFill::operator=(other);
1450  TAttMarker::operator=(other);
1451  fElem = other.fElem;
1452  fElemTop = other.fElemTop;
1453  if (fCoeff) delete [] fCoeff;
1454  fCoeff = 0;
1455  fCsize = other.fCsize;
1456  fNcoeff = other.fNcoeff;
1457  fFactor = other.fFactor;
1458  fTmin = other.fTmin;
1459  fTmax = other.fTmax;
1460  if (fCsize) {
1461  fCoeff = new BtCoef_t[fCsize];
1462  for (Int_t i=0; i<fNcoeff; i++) {
1463  fCoeff[i].cn = other.fCoeff[i].cn;
1464  fCoeff[i].lambda = other.fCoeff[i].lambda;
1465  }
1466  }
1467  return *this;
1468 }
1469 
1470 ////////////////////////////////////////////////////////////////////////////////
1471 /// Addition of other solution.
1472 
1474 {
1475  if (other.GetElement() != fElem) {
1476  Error("operator+=", "Cannot add 2 solutions for different elements");
1477  return *this;
1478  }
1479  Int_t i,j;
1480  BtCoef_t *coeff = fCoeff;
1481  Int_t ncoeff = fNcoeff + other.fNcoeff;
1482  if (ncoeff > fCsize) {
1483  fCsize = ncoeff;
1484  coeff = new BtCoef_t[ncoeff];
1485  for (i=0; i<fNcoeff; i++) {
1486  coeff[i].cn = fCoeff[i].cn;
1487  coeff[i].lambda = fCoeff[i].lambda;
1488  }
1489  delete [] fCoeff;
1490  fCoeff = coeff;
1491  }
1492  ncoeff = fNcoeff;
1493  for (j=0; j<other.fNcoeff; j++) {
1494  for (i=0; i<fNcoeff; i++) {
1495  if (coeff[i].lambda == other.fCoeff[j].lambda) {
1496  coeff[i].cn += fFactor * other.fCoeff[j].cn;
1497  break;
1498  }
1499  }
1500  if (i == fNcoeff) {
1501  coeff[ncoeff].cn = fFactor * other.fCoeff[j].cn;
1502  coeff[ncoeff].lambda = other.fCoeff[j].lambda;
1503  ncoeff++;
1504  }
1505  }
1506  fNcoeff = ncoeff;
1507  return *this;
1508 }
1509 ////////////////////////////////////////////////////////////////////////////////
1510 /// Find concentration of the element at a given time.
1511 
1513 {
1514  Double_t conc = 0.;
1515  for (Int_t i=0; i<fNcoeff; i++)
1516  conc += fCoeff[i].cn * TMath::Exp(-fCoeff[i].lambda * time);
1517  return conc;
1518 }
1519 
1520 ////////////////////////////////////////////////////////////////////////////////
1521 /// Draw the solution of Bateman equation versus time.
1522 
1524 {
1525  if (!gGeoManager) return;
1526  gGeoManager->GetGeomPainter()->DrawBatemanSol(this, option);
1527 }
1528 
1529 ////////////////////////////////////////////////////////////////////////////////
1530 /// Find the solution for the Bateman equations corresponding to the decay
1531 /// chain described by an array ending with element X.
1532 /// A->B->...->X
1533 /// Cn = SUM [Ain * exp(-LMBDi*t)];
1534 /// Cn - concentration Nx/Na
1535 /// n - order of X in chain (A->B->X => n=3)
1536 /// LMBDi - decay constant for element of order i in the chain
1537 /// Ain = LMBD1*...*LMBD(n-1) * br1*...*br(n-1)/(LMBD1-LMBDi)...(LMBDn-LMBDi)
1538 /// bri - branching ratio for decay Ei->Ei+1
1539 
1541 {
1542  fNcoeff = 0;
1543  if (!array || !array->GetEntriesFast()) return;
1544  Int_t n = array->GetEntriesFast();
1545  TGeoDecayChannel *dc = (TGeoDecayChannel*)array->At(n-1);
1546  TGeoElementRN *elem = dc->Daughter();
1547  if (elem != fElem) {
1548  Error("FindSolution", "Last element in the list must be %s\n", fElem->GetName());
1549  return;
1550  }
1551  Int_t i,j;
1552  Int_t order = n+1;
1553  if (!fCoeff) {
1554  fCsize = order;
1555  fCoeff = new BtCoef_t[fCsize];
1556  }
1557  if (fCsize < order) {
1558  delete [] fCoeff;
1559  fCsize = order;
1560  fCoeff = new BtCoef_t[fCsize];
1561  }
1562 
1563  Double_t *lambda = new Double_t[order];
1564  Double_t *br = new Double_t[n];
1565  Double_t halflife;
1566  for (i=0; i<n; i++) {
1567  dc = (TGeoDecayChannel*)array->At(i);
1568  elem = dc->Parent();
1569  br[i] = 0.01 * dc->BranchingRatio();
1570  halflife = elem->HalfLife();
1571  if (halflife==0.) halflife = 1.e-30;
1572  if (elem->Stable()) lambda[i] = 0.;
1573  else lambda[i] = TMath::Log(2.)/halflife;
1574  if (i==n-1) {
1575  elem = dc->Daughter();
1576  halflife = elem->HalfLife();
1577  if (halflife==0.) halflife = 1.e-30;
1578  if (elem->Stable()) lambda[n] = 0.;
1579  else lambda[n] = TMath::Log(2.)/halflife;
1580  }
1581  }
1582  // Check if we have equal lambdas
1583  for (i=0; i<order-1; i++) {
1584  for (j=i+1; j<order; j++) {
1585  if (lambda[j] == lambda[i]) lambda[j] += 0.001*lambda[j];
1586  }
1587  }
1588  Double_t ain;
1589  Double_t pdlambda, plambdabr=1.;
1590  for (j=0; j<n; j++) plambdabr *= lambda[j]*br[j];
1591  for (i=0; i<order; i++) {
1592  pdlambda = 1.;
1593  for (j=0; j<n+1; j++) {
1594  if (j == i) continue;
1595  pdlambda *= lambda[j] - lambda[i];
1596  }
1597  if (pdlambda == 0.) {
1598  Error("FindSolution", "pdlambda=0 !!!");
1599  delete [] lambda;
1600  delete [] br;
1601  return;
1602  }
1603  ain = plambdabr/pdlambda;
1604  fCoeff[i].cn = ain;
1605  fCoeff[i].lambda = lambda[i];
1606  }
1607  fNcoeff = order;
1608  Normalize(fFactor);
1609  delete [] lambda;
1610  delete [] br;
1611 }
1612 
1613 ////////////////////////////////////////////////////////////////////////////////
1614 /// Normalize all coefficients with a given factor.
1615 
1617 {
1618  for (Int_t i=0; i<fNcoeff; i++) fCoeff[i].cn *= factor;
1619 }
1620 
1621 ////////////////////////////////////////////////////////////////////////////////
1622 /// Print concentration evolution.
1623 
1624 void TGeoBatemanSol::Print(Option_t * /*option*/) const
1625 {
1626  TString formula;
1627  formula.Form("N[%s]/N[%s] = ", fElem->GetName(), fElemTop->GetName());
1628  for (Int_t i=0; i<fNcoeff; i++) {
1629  if (i == fNcoeff-1) formula += TString::Format("%g*exp(-%g*t)", fCoeff[i].cn, fCoeff[i].lambda);
1630  else formula += TString::Format("%g*exp(-%g*t) + ", fCoeff[i].cn, fCoeff[i].lambda);
1631  }
1632  printf("%s\n", formula.Data());
1633 }
1634 
static void DecayName(UInt_t decay, TString &name)
Returns decay name.
Int_t GetZ() const
Definition: TGeoElement.h:120
TString fTitle
Definition: TNamed.h:33
TObjArray * GetBranch() const
Definition: TGeoElement.h:352
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
void FillPopulation(TObjArray *population, Double_t precision=0.001, Double_t factor=1.)
Fills the input array with the set of RN elements resulting from the decay of this one...
static TGeoElementRN * ReadElementRN(const char *record, Int_t &ndecays)
Create element from line record.
void Normalize(Double_t factor)
Normalize all coefficients with a given factor.
UnitType setUnitType(UnitType new_type)
Set the currently used unit type (Only ONCE possible)
An array of TObjects.
Definition: TObjArray.h:37
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save primitive for RN elements.
static TGeoDecayChannel * ReadDecay(const char *record)
Create element from line record.
Table of elements.
Definition: TGeoElement.h:369
Double_t fHalfLife
Definition: TGeoElement.h:145
void AddDecay(Int_t decay, Int_t diso, Double_t branchingRatio, Double_t qValue)
Adds a decay mode for this element.
Double_t Log(Double_t x)
Definition: TMath.h:750
Double_t fCoulomb
Definition: TGeoElement.h:51
Double_t fTmin
Definition: TGeoElement.h:293
TLine * line
virtual void Draw(Option_t *option="")
Draw the solution of Bateman equation versus time.
~TGeoBatemanSol()
Destructor.
TGeoElementRN * Daughter() const
Definition: TGeoElement.h:260
const char Option_t
Definition: RtypesCore.h:62
virtual void Print(Option_t *opt=" ") const
Prints decay info.
TObjArray * fListRN
Definition: TGeoElement.h:377
virtual void Print(Option_t *option="") const
Print table of elements.
virtual void Delete(Option_t *option="")
Remove all objects from the array AND delete all heap based objects.
Definition: TObjArray.cxx:355
TGeoElement()
Default constructor.
Definition: TGeoElement.cxx:89
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:687
TGeoElementRN * operator()()
() operator.
static const Int_t gDecayDeltaA[gMaxDecay]
Definition: TGeoElement.cxx:73
TGeoElementRN * fParent
Definition: TGeoElement.h:223
TGeoElementTable * GetElementTable()
Returns material table. Creates it if not existing.
TGeoElementRN * GetElement() const
Definition: TGeoElement.h:310
static const char gElName[gMaxElem][3]
Definition: TGeoElement.cxx:56
void ToUpper()
Change string to upper case.
Definition: TString.cxx:1138
Double_t fFactor
Definition: TGeoElement.h:292
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:172
UnitType unitType()
Access the currently set units type.
void AddRatio(TGeoBatemanSol &ratio)
Adds a proportion ratio to the existing one.
UInt_t Decay() const
Definition: TGeoElement.h:256
Basic string class.
Definition: TString.h:131
Double_t GetA() const
Definition: TGeoElement.h:122
int Int_t
Definition: RtypesCore.h:41
TGeoElementTable & operator=(const TGeoElementTable &)
assignment operator
bool Bool_t
Definition: RtypesCore.h:59
Int_t GetNdecays() const
Get number of decay channels of this element.
void ComputeDerivedQuantities()
Calculate properties for an atomic number.
virtual ~TGeoElementRN()
Destructor.
void AddIsotope(TGeoIsotope *isotope)
Add isotope to the table.
Bool_t CheckDecays() const
Check if all decay chain of the element is OK.
virtual void Print(Option_t *option="") const
This method must be overridden when a class wants to print itself.
Definition: TObject.cxx:550
TObject * At(Int_t idx) const
Definition: TObjArray.h:166
TGeoIsotope * GetIsotope(Int_t i) const
Return i-th isotope in the element.
TGeoIsotope * FindIsotope(const char *name) const
Find existing isotope by name. Not optimized for a big number of isotopes.
Double_t fTH_F
Definition: TGeoElement.h:148
virtual void Print(Option_t *option="") const
Print info about the current decay branch.
virtual void DecayShift(Int_t &dA, Int_t &dZ, Int_t &dI) const
Returns variation in A, Z and Iso after decay.
TObjArray * fIsotopes
Definition: TGeoElement.h:378
TObjArray * fDecays
Definition: TGeoElement.h:155
Double_t * fAbundances
Definition: TGeoElement.h:50
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:694
TGeoElementRN * fElemTop
Definition: TGeoElement.h:289
const TGeoElementRN * fTop
Definition: TGeoElement.h:331
TObjArray * fList
Definition: TGeoElement.h:376
virtual Int_t ENDFCode() const
Definition: TGeoElement.h:178
Marker Attributes class.
Definition: TAttMarker.h:19
TGeoElementRN * Parent() const
Definition: TGeoElement.h:261
static TGeoIsotope * FindIsotope(const char *name)
Find existing isotope by name.
Bool_t HasIsotopes() const
Definition: TGeoElement.h:85
virtual ~TGeoElementTable()
destructor
TGeoBatemanSol & operator=(const TGeoBatemanSol &other)
Assignment.
Fill Area Attributes class.
Definition: TAttFill.h:19
TObjArray * Decays() const
Definition: TGeoElement.h:195
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString...
Definition: TString.cxx:2311
void MakeName(Int_t a, Int_t z, Int_t iso)
Generate a default name for the element.
Double_t fA
Definition: TGeoElement.h:48
Int_t GetN() const
Definition: TGeoElement.h:121
Double_t Neff() const
Returns effective number of nucleons.
Double_t fRadTsai
Definition: TGeoElement.h:52
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save primitive for decays.
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
virtual TObject * FindObject(const char *name) const
Find an object in this collection using its name.
Definition: TObjArray.cxx:414
TObjArray * fBranch
Definition: TGeoElement.h:333
Double_t Concentration(Double_t time) const
Find concentration of the element at a given time.
void AddIsotope(TGeoIsotope *isotope, Double_t relativeAbundance)
Add an isotope for this element. All isotopes have to be isotopes of the same element.
Base class for chemical elements.
Definition: TGeoElement.h:36
static TGeoElementTable * GetElementTable()
Returns pointer to the table.
virtual const char * PrependPathName(const char *dir, TString &name)
Concatenate a directory and a file name.
Definition: TSystem.cxx:1064
static const char gLevName[gMaxLevel]
Definition: TGeoElement.cxx:82
Double_t fTmax
Definition: TGeoElement.h:294
virtual void DrawBatemanSol(TGeoBatemanSol *sol, Option_t *option="")=0
TObject & operator=(const TObject &rhs)
TObject assignment operator.
Definition: TObject.h:268
static constexpr double s
TGeoElementRN()
Default constructor.
virtual ~TGeoElement()
destructor
TGeoElementRN * Down(Int_t ibranch)
Go downwards from current level via ibranch as low in the tree as possible.
BtCoef_t * fCoeff
Definition: TGeoElement.h:295
Double_t HalfLife() const
Definition: TGeoElement.h:186
Bool_t HasRNElements() const
Definition: TGeoElement.h:415
static const Int_t gDecayDeltaZ[gMaxDecay]
Definition: TGeoElement.cxx:78
void SetDaughter(TGeoElementRN *daughter)
Definition: TGeoElement.h:265
Double_t fLimitRatio
Definition: TGeoElement.h:335
TObjArray * fIsotopes
Definition: TGeoElement.h:49
virtual void AddAtAndExpand(TObject *obj, Int_t idx)
Add object at position idx.
Definition: TObjArray.cxx:234
R__EXTERN TSystem * gSystem
Definition: TSystem.h:557
auto * a
Definition: textangle.C:12
static const Int_t gMaxDecay
Definition: TGeoElement.cxx:54
virtual const char * GetName() const
Returns name of decay.
void FindSolution(const TObjArray *array)
Find the solution for the Bateman equations corresponding to the decay chain described by an array en...
virtual TObject * RemoveAt(Int_t idx)
Remove object at index idx.
Definition: TObjArray.cxx:693
const TGeoElementRN * fElem
Definition: TGeoElement.h:332
Double_t fLevel
Definition: TGeoElement.h:143
virtual void Print(Option_t *option="") const
Print info about the element;.
Double_t fTG_F
Definition: TGeoElement.h:149
Double_t fBranchingRatio
Definition: TGeoElement.h:221
void Form(const char *fmt,...)
Formats a string using a printf style format descriptor.
Definition: TString.cxx:2289
unsigned int UInt_t
Definition: RtypesCore.h:42
Int_t GetEntriesFast() const
Definition: TObjArray.h:64
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:880
void AddElementRN(TGeoElementRN *elem)
Add a radionuclide to the table and map it.
Ssiz_t Length() const
Definition: TString.h:405
TGeoDecayChannel & operator=(const TGeoDecayChannel &dc)
Assignment.
static const Int_t gMaxLevel
Definition: TGeoElement.cxx:53
Iterator for decay branches.
Definition: TGeoElement.h:328
TGeoElement * FindElement(const char *name) const
Search an element by symbol or full name Exact matching.
TGeoElementRN * GetElementRN(Int_t ENDFcode) const
Retrieve a radionuclide by ENDF code.
void AddElement(const char *name, const char *title, Int_t z, Double_t a)
Add an element to the table. Obsolete.
static void indent(ostringstream &buf, int indent_level)
TString fName
Definition: TNamed.h:32
constexpr Double_t Na()
Avogadro constant (Avogadro&#39;s Number) in .
Definition: TMath.h:283
TGeoBatemanSol & operator+=(const TGeoBatemanSol &other)
Addition of other solution.
Double_t fDeltaM
Definition: TGeoElement.h:144
static const char * gDecayName[gMaxDecay+1]
Definition: TGeoElement.cxx:68
constexpr Double_t E()
Base of natural log: .
Definition: TMath.h:97
static constexpr double alpha_rcl2
const Bool_t kFALSE
Definition: RtypesCore.h:88
Double_t fNatAbun
Definition: TGeoElement.h:146
Double_t Exp(Double_t x)
Definition: TMath.h:717
static constexpr double fine_structure_const
void SetUsed(Bool_t flag=kTRUE)
Definition: TGeoElement.h:91
Bool_t HasDefaultElements() const
Definition: TGeoElement.h:414
#define ClassImp(name)
Definition: Rtypes.h:365
void ComputeLradTsaiFactor()
Compute Tsai&#39;s Expression for the Radiation Length (Phys Rev. D50 3-1 (1994) page 1254) ...
R__EXTERN TGeoManager * gGeoManager
Definition: TGeoManager.h:600
ElementRNMap_t fElementsRN
Definition: TGeoElement.h:382
Double_t fQvalue
Definition: TGeoElement.h:222
Int_t DecayResult(TGeoDecayChannel *dc) const
Returns ENDF code of decay result.
double Double_t
Definition: RtypesCore.h:55
TGeoElementRN * Next()
Return next element.
Double_t fTH_S
Definition: TGeoElement.h:150
Bool_t CheckTable() const
Checks status of element table.
virtual void Print(Option_t *option="") const
Print this isotope.
Int_t IndexOf(const TObject *obj) const
Definition: TObjArray.cxx:604
static const TString & GetEtcDir()
Get the sysconfig directory in the installation. Static utility function.
Definition: TROOT.cxx:3003
TGeoElementRN * fElem
Definition: TGeoElement.h:288
virtual ~TGeoElemIter()
Destructor.
void SetParent(TGeoElementRN *parent)
Definition: TGeoElement.h:264
void ComputeCoulombFactor()
Compute Coulomb correction factor (Phys Rev. D50 3-1 (1994) page 1254)
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
void SetDefined(Bool_t flag=kTRUE)
Definition: TGeoElement.h:90
Int_t fNisotopes
Definition: TGeoElement.h:47
TGeoIsotope()
Dummy I/O constructor.
Binding & operator=(OUT(*fun)(void))
A decay channel for a radionuclide.
Definition: TGeoElement.h:216
Mother of all ROOT objects.
Definition: TObject.h:37
virtual Double_t GetSpecificActivity() const
Get the activity in Bq of a gram of material made from this element.
you should not use this method at all Int_t Int_t z
Definition: TRolke.cxx:630
Double_t fRatio
Definition: TGeoElement.h:336
virtual void Print(Option_t *option="") const
Print concentration evolution.
TGeoElementRN * Up()
Go upwards from the current location until the next branching, then down.
void ImportElementsRN()
Creates the list of radionuclides.
virtual void Print(Option_t *option="") const
Print this isotope.
Int_t GetIndex() const
Get index of this channel in the list of decays of the parent nuclide.
Class representing a radionuclide.
Definition: TGeoElement.h:138
TGeoElementRN * fDaughter
Definition: TGeoElement.h:224
void BuildDefaultElements()
Creates the default element table.
Int_t GetEntries() const
Return the number of objects in array (i.e.
Definition: TObjArray.cxx:522
TGeoElementTable()
default constructor
static Int_t ENDF(Int_t a, Int_t z, Int_t iso)
Definition: TGeoElement.h:175
TGeoBatemanSol * fRatio
Definition: TGeoElement.h:153
void Add(TObject *obj)
Definition: TObjArray.h:74
Double_t fTG_S
Definition: TGeoElement.h:151
virtual void Fatal(const char *method, const char *msgfmt,...) const
Issue fatal error message.
Definition: TObject.cxx:908
Double_t GetRelativeAbundance(Int_t i) const
Return relative abundance of i-th isotope in this element.
static const Int_t gMaxElem
Definition: TGeoElement.cxx:52
virtual void Print(Option_t *option="") const
Default print for collections, calls Print(option, 1).
Bool_t Stable() const
Definition: TGeoElement.h:194
const Bool_t kTRUE
Definition: RtypesCore.h:87
TVirtualGeoPainter * GetGeomPainter()
Make a default painter if none present. Returns pointer to it.
const Int_t n
Definition: legend1.C:16
Line Attributes class.
Definition: TAttLine.h:18
Double_t BranchingRatio() const
Definition: TGeoElement.h:257
char name[80]
Definition: TGX11.cxx:109
void ResetRatio()
Clears the existing ratio.
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:866
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:48
Double_t fA
Definition: TGeoElement.h:113
const char * Data() const
Definition: TString.h:364
void ExportElementsRN(const char *filename="")
Export radionuclides in a file.
TGeoElemIter & operator=(const TGeoElemIter &iter)
Assignment.