Back to home page

EIC code displayed by LXR

 
 

    


Warning, file /include/root/TGeoElement.h was not indexed or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).

0001 // @(#)root/geom:$Id$
0002 // Author: Andrei Gheata   17/06/04
0003 // Added support for radionuclides: Mihaela Gheata 24/08/2006
0004 /*************************************************************************
0005  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers.               *
0006  * All rights reserved.                                                  *
0007  *                                                                       *
0008  * For the licensing terms see $ROOTSYS/LICENSE.                         *
0009  * For the list of contributors see $ROOTSYS/README/CREDITS.             *
0010  *************************************************************************/
0011 
0012 #ifndef ROOT_TGeoElement
0013 #define ROOT_TGeoElement
0014 
0015 #include "TNamed.h"
0016 
0017 #include "TAttLine.h"
0018 
0019 #include "TAttFill.h"
0020 
0021 #include "TAttMarker.h"
0022 
0023 #include "TObjArray.h"
0024 
0025 #include <map>
0026 
0027 class TGeoElementTable;
0028 class TGeoIsotope;
0029 
0030 ////////////////////////////////////////////////////////////////////////////
0031 //                                                                        //
0032 // TGeoElement - a chemical element                                       //
0033 //                                                                        //
0034 ////////////////////////////////////////////////////////////////////////////
0035 
0036 class TGeoElement : public TNamed {
0037 protected:
0038    enum EGeoElement { kElemUsed = BIT(17), kElemDefined = BIT(18), kElementChecked = BIT(19) };
0039 
0040    Int_t fZ;              // Z of element
0041    Int_t fN;              // Number of nucleons
0042    Int_t fNisotopes;      // Number of isotopes for the element
0043    Double_t fA;           // A of element
0044    TObjArray *fIsotopes;  // List of isotopes
0045    Double_t *fAbundances; //[fNisotopes] Array of relative isotope abundances
0046    Double_t fCoulomb;     // Coulomb correction factor
0047    Double_t fRadTsai;     // Tsai formula for the radiation length
0048 
0049 private:
0050    TGeoElement(const TGeoElement &other) = delete;
0051    TGeoElement &operator=(const TGeoElement &other) = delete;
0052 
0053    // Compute the Coulomb correction factor
0054    void ComputeCoulombFactor();
0055    // Compute the Tsai formula for the radiation length
0056    void ComputeLradTsaiFactor();
0057 
0058 public:
0059    // constructors
0060    TGeoElement();
0061    TGeoElement(const char *name, const char *title, Int_t z, Double_t a);
0062    TGeoElement(const char *name, const char *title, Int_t nisotopes);
0063    TGeoElement(const char *name, const char *title, Int_t z, Int_t n, Double_t a);
0064    // destructor
0065    ~TGeoElement() override;
0066    // methods
0067    virtual Int_t ENDFCode() const { return 0; }
0068    Int_t Z() const { return fZ; }
0069    Int_t N() const { return fN; }
0070    Double_t Neff() const;
0071    Double_t A() const { return fA; }
0072    void AddIsotope(TGeoIsotope *isotope, Double_t relativeAbundance);
0073    Int_t GetNisotopes() const { return fNisotopes; }
0074    TGeoIsotope *GetIsotope(Int_t i) const;
0075    Double_t GetRelativeAbundance(Int_t i) const;
0076    // Calculate properties for an atomic number
0077    void ComputeDerivedQuantities();
0078    // Specific activity (in Bq/gram)
0079    virtual Double_t GetSpecificActivity() const { return 0.; }
0080    Bool_t HasIsotopes() const { return (fNisotopes == 0) ? kFALSE : kTRUE; }
0081    Bool_t IsDefined() const { return TObject::TestBit(kElemDefined); }
0082    virtual Bool_t IsRadioNuclide() const { return kFALSE; }
0083    Bool_t IsUsed() const { return TObject::TestBit(kElemUsed); }
0084    void Print(Option_t *option = "") const override;
0085    void SetDefined(Bool_t flag = kTRUE) { TObject::SetBit(kElemDefined, flag); }
0086    void SetUsed(Bool_t flag = kTRUE) { TObject::SetBit(kElemUsed, flag); }
0087    static TGeoElementTable *GetElementTable();
0088    // Coulomb correction factor
0089    inline Double_t GetfCoulomb() const { return fCoulomb; }
0090    // Tsai formula for the radiation length
0091    inline Double_t GetfRadTsai() const { return fRadTsai; }
0092 
0093    ClassDefOverride(TGeoElement, 3) // base element class
0094 };
0095 
0096 ////////////////////////////////////////////////////////////////////////////
0097 //                                                                        //
0098 // TGeoIsotope - a isotope defined by the atomic number, number of        //
0099 // nucleons and atomic weight (g/mole)                                    //
0100 //                                                                        //
0101 ////////////////////////////////////////////////////////////////////////////
0102 
0103 class TGeoIsotope : public TNamed {
0104 protected:
0105    Int_t fZ;    // atomic number
0106    Int_t fN;    // number of nucleons
0107    Double_t fA; // atomic mass (g/mole)
0108 
0109 public:
0110    TGeoIsotope();
0111    TGeoIsotope(const char *name, Int_t z, Int_t n, Double_t a);
0112    ~TGeoIsotope() override {}
0113 
0114    Int_t GetZ() const { return fZ; }
0115    Int_t GetN() const { return fN; }
0116    Double_t GetA() const { return fA; }
0117    static TGeoIsotope *FindIsotope(const char *name);
0118    void Print(Option_t *option = "") const override;
0119 
0120    ClassDefOverride(TGeoIsotope, 1) // Isotope class defined by Z,N,A
0121 };
0122 
0123 ////////////////////////////////////////////////////////////////////////////
0124 //                                                                        //
0125 // TGeoElementRN - a radionuclide.                                        //
0126 //                                                                        //
0127 ////////////////////////////////////////////////////////////////////////////
0128 
0129 class TGeoDecayChannel;
0130 class TGeoBatemanSol;
0131 
0132 class TGeoElementRN : public TGeoElement {
0133 protected:
0134    Int_t fENDFcode;        // ENDF element code
0135    Int_t fIso;             // Isomer number
0136    Double_t fLevel;        // Isomeric level
0137    Double_t fDeltaM;       // Mass excess
0138    Double_t fHalfLife;     // Half life
0139    Double_t fNatAbun;      // Natural Abundance
0140                            //   char                     fJP[11];   // Spin-parity
0141    Double_t fTH_F;         // Hynalation toxicity
0142    Double_t fTG_F;         // Ingestion toxicity
0143    Double_t fTH_S;         // Hynalation toxicity
0144    Double_t fTG_S;         // Ingestion toxicity
0145    Int_t fStatus;          // Status code
0146    TGeoBatemanSol *fRatio; // Time evolution of proportion by number
0147 
0148    TObjArray *fDecays; // List of decay modes
0149 
0150    void MakeName(Int_t a, Int_t z, Int_t iso);
0151 
0152 private:
0153    TGeoElementRN(const TGeoElementRN &elem) = delete;
0154    TGeoElementRN &operator=(const TGeoElementRN &elem) = delete;
0155 
0156 public:
0157    TGeoElementRN();
0158    TGeoElementRN(Int_t A, Int_t Z, Int_t iso, Double_t level, Double_t deltaM, Double_t halfLife, const char *JP,
0159                  Double_t natAbun, Double_t th_f, Double_t tg_f, Double_t th_s, Double_t tg_s, Int_t status);
0160    ~TGeoElementRN() override;
0161 
0162    void AddDecay(Int_t decay, Int_t diso, Double_t branchingRatio, Double_t qValue);
0163    void AddDecay(TGeoDecayChannel *dc);
0164    void AddRatio(TGeoBatemanSol &ratio);
0165    void ResetRatio();
0166    static Int_t ENDF(Int_t a, Int_t z, Int_t iso) { return 10000 * z + 10 * a + iso; }
0167 
0168    // Getters
0169    Int_t ENDFCode() const override { return fENDFcode; }
0170    Double_t GetSpecificActivity() const override;
0171    Bool_t IsRadioNuclide() const override { return kTRUE; }
0172    Int_t MassNo() const { return (Int_t)fA; }
0173    Int_t AtomicNo() const { return fZ; }
0174    Int_t IsoNo() const { return fIso; }
0175    Double_t Level() const { return fLevel; }
0176    Double_t MassEx() const { return fDeltaM; }
0177    Double_t HalfLife() const { return fHalfLife; }
0178    Double_t NatAbun() const { return fNatAbun; }
0179    const char *PJ() const { return fTitle.Data(); }
0180    Double_t TH_F() const { return fTH_F; }
0181    Double_t TG_F() const { return fTG_F; }
0182    Double_t TH_S() const { return fTH_S; }
0183    Double_t TG_S() const { return fTG_S; }
0184    Double_t Status() const { return fStatus; }
0185    Bool_t Stable() const { return !fDecays; }
0186    TObjArray *Decays() const { return fDecays; }
0187    Int_t GetNdecays() const;
0188    TGeoBatemanSol *Ratio() const { return fRatio; }
0189 
0190    // Utilities
0191    Bool_t CheckDecays() const;
0192    Int_t DecayResult(TGeoDecayChannel *dc) const;
0193    void FillPopulation(TObjArray *population, Double_t precision = 0.001, Double_t factor = 1.);
0194    void Print(Option_t *option = "") const override;
0195    static TGeoElementRN *ReadElementRN(const char *record, Int_t &ndecays);
0196    void SavePrimitive(std::ostream &out, Option_t *option = "") override;
0197 
0198    ClassDefOverride(TGeoElementRN, 2) // radionuclides class
0199 };
0200 
0201 ////////////////////////////////////////////////////////////////////////////
0202 //                                                                        //
0203 // TGeoDecayChannel - decay channel utility class.                        //
0204 //                                                                        //
0205 ////////////////////////////////////////////////////////////////////////////
0206 
0207 class TGeoDecayChannel : public TObject {
0208 private:
0209    UInt_t fDecay;            // Decay mode
0210    Int_t fDiso;              // Delta isomeric number
0211    Double_t fBranchingRatio; // Branching Ratio
0212    Double_t fQvalue;         // Qvalue in GeV
0213    TGeoElementRN *fParent;   // Parent element
0214    TGeoElementRN *fDaughter; // Daughter element
0215 public:
0216    enum ENuclearDecayMode {
0217       kBitMask32 = 0xffffffff,
0218       k2BetaMinus = BIT(0),
0219       kBetaMinus = BIT(1),
0220       kNeutronEm = BIT(2),
0221       kProtonEm = BIT(3),
0222       kAlpha = BIT(4),
0223       kECF = BIT(5),
0224       kElecCapt = BIT(6),
0225       kIsoTrans = BIT(7),
0226       kI = BIT(8),
0227       kSpontFiss = BIT(9),
0228       k2P = BIT(10),
0229       k2N = BIT(11),
0230       k2A = BIT(12),
0231       kCarbon12 = BIT(13),
0232       kCarbon14 = BIT(14)
0233    };
0234    TGeoDecayChannel() : fDecay(0), fDiso(0), fBranchingRatio(0), fQvalue(0), fParent(nullptr), fDaughter(nullptr) {}
0235    TGeoDecayChannel(Int_t decay, Int_t diso, Double_t branchingRatio, Double_t qValue)
0236       : fDecay(decay),
0237         fDiso(diso),
0238         fBranchingRatio(branchingRatio),
0239         fQvalue(qValue),
0240         fParent(nullptr),
0241         fDaughter(nullptr)
0242    {
0243    }
0244    TGeoDecayChannel(const TGeoDecayChannel &dc)
0245       : TObject(dc),
0246         fDecay(dc.fDecay),
0247         fDiso(dc.fDiso),
0248         fBranchingRatio(dc.fBranchingRatio),
0249         fQvalue(dc.fQvalue),
0250         fParent(dc.fParent),
0251         fDaughter(dc.fDaughter)
0252    {
0253    }
0254    ~TGeoDecayChannel() override {}
0255 
0256    TGeoDecayChannel &operator=(const TGeoDecayChannel &dc);
0257 
0258    // Getters
0259    Int_t GetIndex() const;
0260    const char *GetName() const override;
0261    UInt_t Decay() const { return fDecay; }
0262    Double_t BranchingRatio() const { return fBranchingRatio; }
0263    Double_t Qvalue() const { return fQvalue; }
0264    Int_t DeltaIso() const { return fDiso; }
0265    TGeoElementRN *Daughter() const { return fDaughter; }
0266    TGeoElementRN *Parent() const { return fParent; }
0267    static void DecayName(UInt_t decay, TString &name);
0268    // Setters
0269    void SetParent(TGeoElementRN *parent) { fParent = parent; }
0270    void SetDaughter(TGeoElementRN *daughter) { fDaughter = daughter; }
0271    // Services
0272    void Print(Option_t *opt = " ") const override;
0273    static TGeoDecayChannel *ReadDecay(const char *record);
0274    void SavePrimitive(std::ostream &out, Option_t *option = "") override;
0275    virtual void DecayShift(Int_t &dA, Int_t &dZ, Int_t &dI) const;
0276 
0277    ClassDefOverride(TGeoDecayChannel, 1) // Decay channel for Elements
0278 };
0279 
0280 ////////////////////////////////////////////////////////////////////////////////
0281 //                                                                            //
0282 // TGeoBatemanSol -Class representing the Bateman solution for a decay branch //
0283 //                                                                            //
0284 ////////////////////////////////////////////////////////////////////////////////
0285 
0286 class TGeoBatemanSol : public TObject, public TAttLine, public TAttFill, public TAttMarker {
0287 private:
0288    typedef struct {
0289       Double_t cn;     // Concentration for element 'i': Ni/Ntop
0290       Double_t lambda; // Decay coef. for element 'i'
0291    } BtCoef_t;
0292    TGeoElementRN *fElem;    // Referred RN element
0293    TGeoElementRN *fElemTop; // Top RN element
0294    Int_t fCsize;            // Size of the array of coefficients
0295    Int_t fNcoeff;           // Number of coefficients
0296    Double_t fFactor;        // Constant factor that applies to all coefficients
0297    Double_t fTmin;          // Minimum value of the time interval
0298    Double_t fTmax;          // Maximum value of the time interval
0299    BtCoef_t *fCoeff;        //[fNcoeff] Array of coefficients
0300 public:
0301    TGeoBatemanSol()
0302       : TObject(),
0303         TAttLine(),
0304         TAttFill(),
0305         TAttMarker(),
0306         fElem(nullptr),
0307         fElemTop(nullptr),
0308         fCsize(0),
0309         fNcoeff(0),
0310         fFactor(1.),
0311         fTmin(0.),
0312         fTmax(0),
0313         fCoeff(nullptr)
0314    {
0315    }
0316    TGeoBatemanSol(TGeoElementRN *elem);
0317    TGeoBatemanSol(const TObjArray *chain);
0318    TGeoBatemanSol(const TGeoBatemanSol &other);
0319    ~TGeoBatemanSol() override;
0320 
0321    TGeoBatemanSol &operator=(const TGeoBatemanSol &other);
0322    TGeoBatemanSol &operator+=(const TGeoBatemanSol &other);
0323 
0324    Double_t Concentration(Double_t time) const;
0325    void Draw(Option_t *option = "") override;
0326    void GetCoeff(Int_t i, Double_t &cn, Double_t &lambda) const
0327    {
0328       cn = fCoeff[i].cn;
0329       lambda = fCoeff[i].lambda;
0330    }
0331    void GetRange(Double_t &tmin, Double_t &tmax) const
0332    {
0333       tmin = fTmin;
0334       tmax = fTmax;
0335    }
0336    TGeoElementRN *GetElement() const { return fElem; }
0337    TGeoElementRN *GetTopElement() const { return fElemTop; }
0338    Int_t GetNcoeff() const { return fNcoeff; }
0339    void Print(Option_t *option = "") const override;
0340    void SetRange(Double_t tmin = 0., Double_t tmax = 0.)
0341    {
0342       fTmin = tmin;
0343       fTmax = tmax;
0344    }
0345    void SetFactor(Double_t factor) { fFactor = factor; }
0346    void FindSolution(const TObjArray *array);
0347    void Normalize(Double_t factor);
0348 
0349    ClassDefOverride(TGeoBatemanSol, 1) // Solution for the Bateman equation
0350 };
0351 
0352 ////////////////////////////////////////////////////////////////////////////
0353 //                                                                        //
0354 // TGeoElemIter - iterator for decay chains.                              //
0355 //                                                                        //
0356 ////////////////////////////////////////////////////////////////////////////
0357 
0358 class TGeoElemIter {
0359 private:
0360    const TGeoElementRN *fTop;  // Top element of the iteration
0361    const TGeoElementRN *fElem; // Current element
0362    TObjArray *fBranch;         // Current branch
0363    Int_t fLevel;               // Current level
0364    Double_t fLimitRatio;       // Minimum cumulative branching ratio
0365    Double_t fRatio;            // Current ratio
0366 
0367 protected:
0368    TGeoElemIter() : fTop(nullptr), fElem(nullptr), fBranch(nullptr), fLevel(0), fLimitRatio(0), fRatio(0) {}
0369    TGeoElementRN *Down(Int_t ibranch);
0370    TGeoElementRN *Up();
0371 
0372 public:
0373    TGeoElemIter(TGeoElementRN *top, Double_t limit = 1.e-4);
0374    TGeoElemIter(const TGeoElemIter &iter);
0375    virtual ~TGeoElemIter();
0376 
0377    TGeoElemIter &operator=(const TGeoElemIter &iter);
0378    TGeoElementRN *operator()();
0379    TGeoElementRN *Next();
0380 
0381    TObjArray *GetBranch() const { return fBranch; }
0382    const TGeoElementRN *GetTop() const { return fTop; }
0383    const TGeoElementRN *GetElement() const { return fElem; }
0384    Int_t GetLevel() const { return fLevel; }
0385    Double_t GetRatio() const { return fRatio; }
0386    virtual void Print(Option_t *option = "") const;
0387    void SetLimitRatio(Double_t limit) { fLimitRatio = limit; }
0388 
0389    ClassDef(TGeoElemIter, 0) // Iterator for radionuclide chains.
0390 };
0391 
0392 ////////////////////////////////////////////////////////////////////////////
0393 //                                                                        //
0394 // TGeoElementTable - table of elements                                   //
0395 //                                                                        //
0396 ////////////////////////////////////////////////////////////////////////////
0397 
0398 class TGeoElementTable : public TObject {
0399 private:
0400    // data members
0401    Int_t fNelements;     // number of elements
0402    Int_t fNelementsRN;   // number of RN elements
0403    Int_t fNisotopes;     // number of isotopes
0404    TObjArray *fList;     // list of elements
0405    TObjArray *fListRN;   // list of RN elements
0406    TObjArray *fIsotopes; // list of user-defined isotopes
0407    // Map of radionuclides
0408    typedef std::map<Int_t, TGeoElementRN *> ElementRNMap_t;
0409    typedef ElementRNMap_t::iterator ElementRNMapIt_t;
0410    ElementRNMap_t fElementsRN; //! map of RN elements with ENDF key
0411 
0412 protected:
0413    TGeoElementTable(const TGeoElementTable &);
0414    TGeoElementTable &operator=(const TGeoElementTable &);
0415 
0416 public:
0417    // constructors
0418    TGeoElementTable();
0419    TGeoElementTable(Int_t nelements);
0420    // destructor
0421    ~TGeoElementTable() override;
0422    // methods
0423 
0424    enum EGeoETStatus { kETDefaultElements = BIT(14), kETRNElements = BIT(15) };
0425    void AddElement(const char *name, const char *title, Int_t z, Double_t a);
0426    void AddElement(const char *name, const char *title, Int_t z, Int_t n, Double_t a);
0427    void AddElement(TGeoElement *elem);
0428    void AddElementRN(TGeoElementRN *elem);
0429    void AddIsotope(TGeoIsotope *isotope);
0430    void BuildDefaultElements();
0431    void ImportElementsRN();
0432    Bool_t CheckTable() const;
0433    TGeoElement *FindElement(const char *name) const;
0434    TGeoIsotope *FindIsotope(const char *name) const;
0435    TGeoElement *GetElement(Int_t z) { return (TGeoElement *)fList->At(z); }
0436    TGeoElementRN *GetElementRN(Int_t ENDFcode) const;
0437    TGeoElementRN *GetElementRN(Int_t a, Int_t z, Int_t iso = 0) const;
0438    TObjArray *GetElementsRN() const { return fListRN; }
0439    Bool_t HasDefaultElements() const { return TObject::TestBit(kETDefaultElements); }
0440    Bool_t HasRNElements() const { return TObject::TestBit(kETRNElements); }
0441 
0442    Int_t GetNelements() const { return fNelements; }
0443    Int_t GetNelementsRN() const { return fNelementsRN; }
0444    void ExportElementsRN(const char *filename = "");
0445    void Print(Option_t *option = "") const override;
0446 
0447    ClassDefOverride(TGeoElementTable, 4) // table of elements
0448 };
0449 
0450 #endif