Logo ROOT   6.30.04
Reference Guide
 All Namespaces Files Pages
TGeoElement.h
Go to the documentation of this file.
1 // @(#)root/geom:$Id$
2 // Author: Andrei Gheata 17/06/04
3 // Added support for radionuclides: Mihaela Gheata 24/08/2006
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 #ifndef ROOT_TGeoElement
13 #define ROOT_TGeoElement
14 
15 #include "TNamed.h"
16 
17 #include "TAttLine.h"
18 
19 #include "TAttFill.h"
20 
21 #include "TAttMarker.h"
22 
23 #include "TObjArray.h"
24 
25 #include <map>
26 
27 class TGeoElementTable;
28 class TGeoIsotope;
29 
30 ////////////////////////////////////////////////////////////////////////////
31 // //
32 // TGeoElement - a chemical element //
33 // //
34 ////////////////////////////////////////////////////////////////////////////
35 
36 class TGeoElement : public TNamed
37 {
38 protected:
39  enum EGeoElement {
40  kElemUsed = BIT(17),
41  kElemDefined = BIT(18),
42  kElementChecked = BIT(19)
43  };
44 
45  Int_t fZ; // Z of element
46  Int_t fN; // Number of nucleons
47  Int_t fNisotopes; // Number of isotopes for the element
48  Double_t fA; // A of element
49  TObjArray *fIsotopes; // List of isotopes
50  Double_t *fAbundances; //[fNisotopes] Array of relative isotope abundances
51  Double_t fCoulomb; // Coulomb correction factor
52  Double_t fRadTsai; // Tsai formula for the radiation length
53 
54 private:
55  TGeoElement(const TGeoElement &other);
56  TGeoElement &operator=(const TGeoElement &other);
57 
58  // Compute the Coulomb correction factor
59  void ComputeCoulombFactor();
60  // Compute the Tsai formula for the radiation length
61  void ComputeLradTsaiFactor();
62 
63 public:
64  // constructors
65  TGeoElement();
66  TGeoElement(const char *name, const char *title, Int_t z, Double_t a);
67  TGeoElement(const char *name, const char *title, Int_t nisotopes);
68  TGeoElement(const char *name, const char *title, Int_t z, Int_t n, Double_t a);
69  // destructor
70  virtual ~TGeoElement() {;}
71  // methods
72  virtual Int_t ENDFCode() const { return 0;}
73  Int_t Z() const {return fZ;}
74  Int_t N() const {return fN;}
75  Double_t Neff() const;
76  Double_t A() const {return fA;}
77  void AddIsotope(TGeoIsotope *isotope, Double_t relativeAbundance);
78  Int_t GetNisotopes() const {return fNisotopes;}
79  TGeoIsotope *GetIsotope(Int_t i) const;
80  Double_t GetRelativeAbundance(Int_t i) const;
81  // Calculate properties for an atomic number
82  void ComputeDerivedQuantities();
83  // Specific activity (in Bq/gram)
84  virtual Double_t GetSpecificActivity() const {return 0.;}
85  Bool_t HasIsotopes() const {return (fNisotopes==0)?kFALSE:kTRUE;}
86  Bool_t IsDefined() const {return TObject::TestBit(kElemDefined);}
87  virtual Bool_t IsRadioNuclide() const {return kFALSE;}
88  Bool_t IsUsed() const {return TObject::TestBit(kElemUsed);}
89  virtual void Print(Option_t *option = "") const;
90  void SetDefined(Bool_t flag=kTRUE) {TObject::SetBit(kElemDefined,flag);}
91  void SetUsed(Bool_t flag=kTRUE) {TObject::SetBit(kElemUsed,flag);}
92  static TGeoElementTable *GetElementTable();
93  // Coulomb correction factor
94  inline Double_t GetfCoulomb() const {return fCoulomb;}
95  // Tsai formula for the radiation length
96  inline Double_t GetfRadTsai() const {return fRadTsai;}
97 
98  ClassDef(TGeoElement, 3) // base element class
99 };
100 
101 ////////////////////////////////////////////////////////////////////////////
102 // //
103 // TGeoIsotope - a isotope defined by the atomic number, number of //
104 // nucleons and atomic weight (g/mole) //
105 // //
106 ////////////////////////////////////////////////////////////////////////////
107 
108 class TGeoIsotope : public TNamed
109 {
110 protected:
111  Int_t fZ; // atomic number
112  Int_t fN; // number of nucleons
113  Double_t fA; // atomic mass (g/mole)
114 
115 public:
116  TGeoIsotope();
117  TGeoIsotope(const char *name, Int_t z, Int_t n, Double_t a);
118  virtual ~TGeoIsotope() {}
119 
120  Int_t GetZ() const {return fZ;}
121  Int_t GetN() const {return fN;}
122  Double_t GetA() const {return fA;}
123  static TGeoIsotope *FindIsotope(const char *name);
124  virtual void Print(Option_t *option = "") const;
125 
126  ClassDef(TGeoIsotope, 1) // Isotope class defined by Z,N,A
127 };
128 
129 ////////////////////////////////////////////////////////////////////////////
130 // //
131 // TGeoElementRN - a radionuclide. //
132 // //
133 ////////////////////////////////////////////////////////////////////////////
134 
135 class TGeoDecayChannel;
136 class TGeoBatemanSol;
137 
138 class TGeoElementRN : public TGeoElement
139 {
140 protected:
141  Int_t fENDFcode; // ENDF element code
142  Int_t fIso; // Isomer number
143  Double_t fLevel; // Isomeric level
144  Double_t fDeltaM; // Mass excess
145  Double_t fHalfLife; // Half life
146  Double_t fNatAbun; // Natural Abundance
147 // char fJP[11]; // Spin-parity
148  Double_t fTH_F; // Hynalation toxicity
149  Double_t fTG_F; // Ingestion toxicity
150  Double_t fTH_S; // Hynalation toxicity
151  Double_t fTG_S; // Ingestion toxicity
152  Int_t fStatus; // Status code
153  TGeoBatemanSol *fRatio; // Time evolution of proportion by number
154 
155  TObjArray *fDecays; // List of decay modes
156 
157  void MakeName(Int_t a, Int_t z, Int_t iso);
158 
159 private:
160  TGeoElementRN(const TGeoElementRN& elem);
161  TGeoElementRN& operator=(const TGeoElementRN& elem);
162 
163 public:
164  TGeoElementRN();
165  TGeoElementRN(Int_t A, Int_t Z, Int_t iso, Double_t level,
166  Double_t deltaM, Double_t halfLife, const char* JP,
167  Double_t natAbun, Double_t th_f, Double_t tg_f, Double_t th_s,
168  Double_t tg_s, Int_t status);
169  virtual ~TGeoElementRN();
170 
171  void AddDecay(Int_t decay, Int_t diso, Double_t branchingRatio, Double_t qValue);
172  void AddDecay(TGeoDecayChannel *dc);
173  void AddRatio(TGeoBatemanSol &ratio);
174  void ResetRatio();
175  static Int_t ENDF(Int_t a, Int_t z, Int_t iso) {return 10000*z+10*a+iso;}
176 
177  // Getters
178  virtual Int_t ENDFCode() const {return fENDFcode;}
179  virtual Double_t GetSpecificActivity() const;
180  virtual Bool_t IsRadioNuclide() const {return kTRUE;}
181  Int_t MassNo() const {return (Int_t)fA;}
182  Int_t AtomicNo() const {return fZ;}
183  Int_t IsoNo() const {return fIso;}
184  Double_t Level() const {return fLevel;}
185  Double_t MassEx() const {return fDeltaM;}
186  Double_t HalfLife() const {return fHalfLife;}
187  Double_t NatAbun() const {return fNatAbun;}
188  const char* PJ() const {return fTitle.Data();}
189  Double_t TH_F() const {return fTH_F;}
190  Double_t TG_F() const {return fTG_F;}
191  Double_t TH_S() const {return fTH_S;}
192  Double_t TG_S() const {return fTG_S;}
193  Double_t Status() const {return fStatus;}
194  Bool_t Stable() const {return !fDecays;}
195  TObjArray *Decays() const {return fDecays;}
196  Int_t GetNdecays() const;
197  TGeoBatemanSol *Ratio() const {return fRatio;}
198 
199  // Utilities
200  Bool_t CheckDecays() const;
201  Int_t DecayResult(TGeoDecayChannel *dc) const;
202  void FillPopulation(TObjArray *population, Double_t precision=0.001, Double_t factor=1.);
203  virtual void Print(Option_t *option = "") const;
204  static TGeoElementRN *ReadElementRN(const char *record, Int_t &ndecays);
205  virtual void SavePrimitive(std::ostream &out, Option_t *option = "");
206 
207  ClassDef(TGeoElementRN, 2) // radionuclides class
208 };
209 
210 ////////////////////////////////////////////////////////////////////////////
211 // //
212 // TGeoDecayChannel - decay channel utility class. //
213 // //
214 ////////////////////////////////////////////////////////////////////////////
215 
216 class TGeoDecayChannel : public TObject
217 {
218 private:
219  UInt_t fDecay; // Decay mode
220  Int_t fDiso; // Delta isomeric number
221  Double_t fBranchingRatio; // Branching Ratio
222  Double_t fQvalue; // Qvalue in GeV
223  TGeoElementRN *fParent; // Parent element
224  TGeoElementRN *fDaughter; // Daughter element
225 public:
226  enum ENuclearDecayMode {
227  kBitMask32 = 0xffffffff,
228  k2BetaMinus = BIT(0),
229  kBetaMinus = BIT(1),
230  kNeutronEm = BIT(2),
231  kProtonEm = BIT(3),
232  kAlpha = BIT(4),
233  kECF = BIT(5),
234  kElecCapt = BIT(6),
235  kIsoTrans = BIT(7),
236  kI = BIT(8),
237  kSpontFiss = BIT(9),
238  k2P = BIT(10),
239  k2N = BIT(11),
240  k2A = BIT(12),
241  kCarbon12 = BIT(13),
242  kCarbon14 = BIT(14)
243  };
244  TGeoDecayChannel() : fDecay(0), fDiso(0), fBranchingRatio(0), fQvalue(0), fParent(0), fDaughter(0) {}
245  TGeoDecayChannel(Int_t decay, Int_t diso, Double_t branchingRatio, Double_t qValue)
246  : fDecay(decay), fDiso(diso), fBranchingRatio(branchingRatio), fQvalue(qValue), fParent(0), fDaughter(0) {}
247  TGeoDecayChannel(const TGeoDecayChannel &dc) : TObject(dc),fDecay(dc.fDecay),fDiso(dc.fDiso),fBranchingRatio(dc.fBranchingRatio),
248  fQvalue(dc.fQvalue),fParent(dc.fParent),fDaughter(dc.fDaughter) {}
249  virtual ~TGeoDecayChannel() {}
250 
251  TGeoDecayChannel& operator=(const TGeoDecayChannel& dc);
252 
253  // Getters
254  Int_t GetIndex() const;
255  virtual const char *GetName() const;
256  UInt_t Decay() const {return fDecay;}
257  Double_t BranchingRatio() const {return fBranchingRatio;}
258  Double_t Qvalue() const {return fQvalue;}
259  Int_t DeltaIso() const {return fDiso;}
260  TGeoElementRN *Daughter() const {return fDaughter;}
261  TGeoElementRN *Parent() const {return fParent;}
262  static void DecayName(UInt_t decay, TString &name);
263  // Setters
264  void SetParent(TGeoElementRN *parent) {fParent = parent;}
265  void SetDaughter(TGeoElementRN *daughter) {fDaughter = daughter;}
266  // Services
267  virtual void Print(Option_t *opt = " ") const;
268  static TGeoDecayChannel *ReadDecay(const char *record);
269  virtual void SavePrimitive(std::ostream &out, Option_t *option = "");
270  virtual void DecayShift(Int_t &dA, Int_t &dZ, Int_t &dI) const ;
271 
272  ClassDef(TGeoDecayChannel,1) // Decay channel for Elements
273 };
274 
275 ////////////////////////////////////////////////////////////////////////////////
276 // //
277 // TGeoBatemanSol -Class representing the Bateman solution for a decay branch //
278 // //
279 ////////////////////////////////////////////////////////////////////////////////
280 
281 class TGeoBatemanSol : public TObject, public TAttLine, public TAttFill, public TAttMarker
282 {
283 private:
284  typedef struct {
285  Double_t cn; // Concentration for element 'i': Ni/Ntop
286  Double_t lambda; // Decay coef. for element 'i'
287  } BtCoef_t;
288  TGeoElementRN *fElem; // Referred RN element
289  TGeoElementRN *fElemTop; // Top RN element
290  Int_t fCsize; // Size of the array of coefficients
291  Int_t fNcoeff; // Number of coefficients
292  Double_t fFactor; // Constant factor that applies to all coefficients
293  Double_t fTmin; // Minimum value of the time interval
294  Double_t fTmax; // Maximum value of the time interval
295  BtCoef_t *fCoeff; //[fNcoeff] Array of coefficients
296 public:
297  TGeoBatemanSol() : TObject(), TAttLine(), TAttFill(), TAttMarker(), fElem(NULL), fElemTop(NULL), fCsize(0), fNcoeff(0), fFactor(1.), fTmin(0.), fTmax(0), fCoeff(NULL) {}
298  TGeoBatemanSol(TGeoElementRN *elem);
299  TGeoBatemanSol(const TObjArray *chain);
300  TGeoBatemanSol(const TGeoBatemanSol& other);
301  ~TGeoBatemanSol();
302 
303  TGeoBatemanSol& operator=(const TGeoBatemanSol& other);
304  TGeoBatemanSol& operator+=(const TGeoBatemanSol& other);
305 
306  Double_t Concentration(Double_t time) const;
307  virtual void Draw(Option_t *option="");
308  void GetCoeff(Int_t i, Double_t &cn, Double_t &lambda) const {cn=fCoeff[i].cn; lambda=fCoeff[i].lambda;}
309  void GetRange(Double_t &tmin, Double_t &tmax) const {tmin=fTmin; tmax=fTmax;}
310  TGeoElementRN *GetElement() const {return fElem;}
311  TGeoElementRN *GetTopElement() const {return fElemTop;}
312  Int_t GetNcoeff() const {return fNcoeff;}
313  virtual void Print(Option_t *option = "") const;
314  void SetRange(Double_t tmin=0., Double_t tmax=0.) {fTmin=tmin; fTmax=tmax;}
315  void SetFactor(Double_t factor) {fFactor = factor;}
316  void FindSolution(const TObjArray *array);
317  void Normalize(Double_t factor);
318 
319  ClassDef(TGeoBatemanSol,1) // Solution for the Bateman equation
320 };
321 
322 ////////////////////////////////////////////////////////////////////////////
323 // //
324 // TGeoElemIter - iterator for decay chains. //
325 // //
326 ////////////////////////////////////////////////////////////////////////////
327 
328 class TGeoElemIter
329 {
330 private:
331  const TGeoElementRN *fTop; // Top element of the iteration
332  const TGeoElementRN *fElem; // Current element
333  TObjArray *fBranch; // Current branch
334  Int_t fLevel; // Current level
335  Double_t fLimitRatio; // Minimum cumulative branching ratio
336  Double_t fRatio; // Current ratio
337 
338 protected:
339  TGeoElemIter() : fTop(0), fElem(0), fBranch(0), fLevel(0), fLimitRatio(0), fRatio(0) {}
340  TGeoElementRN *Down(Int_t ibranch);
341  TGeoElementRN *Up();
342 
343 public:
344  TGeoElemIter(TGeoElementRN *top, Double_t limit=1.e-4);
345  TGeoElemIter(const TGeoElemIter &iter);
346  virtual ~TGeoElemIter();
347 
348  TGeoElemIter &operator=(const TGeoElemIter &iter);
349  TGeoElementRN *operator()();
350  TGeoElementRN *Next();
351 
352  TObjArray *GetBranch() const {return fBranch;}
353  const TGeoElementRN *GetTop() const {return fTop;}
354  const TGeoElementRN *GetElement() const {return fElem;}
355  Int_t GetLevel() const {return fLevel;}
356  Double_t GetRatio() const {return fRatio;}
357  virtual void Print(Option_t *option="") const;
358  void SetLimitRatio(Double_t limit) {fLimitRatio = limit;}
359 
360  ClassDef(TGeoElemIter,0) // Iterator for radionuclide chains.
361 };
362 
363 ////////////////////////////////////////////////////////////////////////////
364 // //
365 // TGeoElementTable - table of elements //
366 // //
367 ////////////////////////////////////////////////////////////////////////////
368 
369 class TGeoElementTable : public TObject
370 {
371 private:
372 // data members
373  Int_t fNelements; // number of elements
374  Int_t fNelementsRN; // number of RN elements
375  Int_t fNisotopes; // number of isotopes
376  TObjArray *fList; // list of elements
377  TObjArray *fListRN; // list of RN elements
378  TObjArray *fIsotopes; // list of user-defined isotopes
379  // Map of radionuclides
380  typedef std::map<Int_t, TGeoElementRN *> ElementRNMap_t;
381  typedef ElementRNMap_t::iterator ElementRNMapIt_t;
382  ElementRNMap_t fElementsRN; //! map of RN elements with ENDF key
383 
384 protected:
385  TGeoElementTable(const TGeoElementTable&);
386  TGeoElementTable& operator=(const TGeoElementTable&);
387 
388 public:
389  // constructors
390  TGeoElementTable();
391  TGeoElementTable(Int_t nelements);
392  // destructor
393  virtual ~TGeoElementTable();
394  // methods
395 
396  enum EGeoETStatus {
397  kETDefaultElements = BIT(14),
398  kETRNElements = BIT(15)
399  };
400  void AddElement(const char *name, const char *title, Int_t z, Double_t a);
401  void AddElement(const char *name, const char *title, Int_t z, Int_t n, Double_t a);
402  void AddElement(TGeoElement *elem);
403  void AddElementRN(TGeoElementRN *elem);
404  void AddIsotope(TGeoIsotope *isotope);
405  void BuildDefaultElements();
406  void ImportElementsRN();
407  Bool_t CheckTable() const;
408  TGeoElement *FindElement(const char *name) const;
409  TGeoIsotope *FindIsotope(const char *name) const;
410  TGeoElement *GetElement(Int_t z) {return (TGeoElement*)fList->At(z);}
411  TGeoElementRN *GetElementRN(Int_t ENDFcode) const;
412  TGeoElementRN *GetElementRN(Int_t a, Int_t z, Int_t iso=0) const;
413  TObjArray *GetElementsRN() const {return fListRN;}
414  Bool_t HasDefaultElements() const {return TObject::TestBit(kETDefaultElements);}
415  Bool_t HasRNElements() const {return TObject::TestBit(kETRNElements);}
416 
417  Int_t GetNelements() const {return fNelements;}
418  Int_t GetNelementsRN() const {return fNelementsRN;}
419  void ExportElementsRN(const char *filename="");
420  virtual void Print(Option_t *option = "") const;
421 
422  ClassDef(TGeoElementTable,4) // table of elements
423 };
424 
425 #endif
426