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