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