Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TGeoMaterial.cxx
Go to the documentation of this file.
1// @(#)root/geom:$Id$
2// Author: Andrei Gheata 25/10/01
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 TGeoMaterial
13\ingroup Materials_classes
14
15Base class describing materials.
16
17## Important note about units
18Since **v6-17-02** the geometry package adopted a system of units, upon the request to support
19an in-memory material representation consistent with the one in Geant4. The adoption was done
20gradually and starting with **v6-19-02** (back-ported to **v6-18-02**) the package supports changing
21the default units to either ROOT (CGS) or Geant4 ones. In the same version the Geant4 units were
22set to be the default ones, changing the previous behavior and making material properties such
23as radiation and interaction lengths having in memory values an order of magnitude lower. This behavior
24affected versions up to **v6-25-01**, after which the default units were restored to be the ROOT ones.
25
26For users needing to restore the CGS behavior for material properties, the following sequence needs
27to be called before creating the TGeoManager instance:
28 * From **v6-18-02** to **v6-22-06**:
29```
30 TGeoUnit::setUnitType(TGeoUnit::kTGeoUnits);
31```
32
33 * From **v6-22-08** to **v6-25-01**:
34```
35 TGeoManager::LockDefaultUnits(false);
36 TGeoManager::SetDefaultUnits(kRootUnits);
37 TGeoManager::LockDefaultUnits(true);
38```
39*/
40
41#include <iostream>
42#include <limits>
43#include "TMath.h"
44#include "TObjArray.h"
45#include "TGeoElement.h"
46#include "TGeoManager.h"
47#include "TGeoExtension.h"
48#include "TGeoMaterial.h"
51#include "TGDMLMatrix.h"
52
53// statics and globals
54
55
56////////////////////////////////////////////////////////////////////////////////
57/// Default constructor.
58
60 : TNamed(),
61 TAttFill(),
62 fIndex(0),
63 fA(0.),
64 fZ(0.),
65 fDensity(0.),
66 fRadLen(0.),
67 fIntLen(0.),
68 fTemperature(0.),
69 fPressure(0.),
70 fState(kMatStateUndefined),
71 fShader(nullptr),
72 fCerenkov(nullptr),
73 fElement(nullptr),
74 fUserExtension(nullptr),
75 fFWExtension(nullptr)
76{
77 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
79 fIndex = -1;
83}
84
85////////////////////////////////////////////////////////////////////////////////
86/// Constructor.
87///
88/// \param name material name.
89
91 : TNamed(name, ""),
92 TAttFill(),
93 fIndex(0),
94 fA(0.),
95 fZ(0.),
96 fDensity(0.),
97 fRadLen(0.),
98 fIntLen(0.),
99 fTemperature(0.),
100 fPressure(0.),
101 fState(kMatStateUndefined),
102 fShader(nullptr),
103 fCerenkov(nullptr),
104 fElement(nullptr),
105 fUserExtension(nullptr),
106 fFWExtension(nullptr)
107{
108 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
109 fName = fName.Strip();
111 fIndex = -1;
115
116 if (!gGeoManager) {
117 gGeoManager = new TGeoManager("Geometry", "default geometry");
118 }
120}
121
122////////////////////////////////////////////////////////////////////////////////
123/// Constructor.
124///
125/// \param name material name.
126/// \param a atomic mass.
127/// \param z atomic number.
128/// \param rho material density in g/cm3.
129/// \param radlen
130/// \param intlen
131
133 : TNamed(name, ""),
134 TAttFill(),
135 fIndex(0),
136 fA(a),
137 fZ(z),
138 fDensity(rho),
139 fRadLen(0.),
140 fIntLen(0.),
141 fTemperature(0.),
142 fPressure(0.),
143 fState(kMatStateUndefined),
144 fShader(nullptr),
145 fCerenkov(nullptr),
146 fElement(nullptr),
147 fUserExtension(nullptr),
148 fFWExtension(nullptr)
149{
150 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
151 fName = fName.Strip();
153 fIndex = -1;
154 fA = a;
155 fZ = z;
156 fDensity = rho;
161 if (!gGeoManager) {
162 gGeoManager = new TGeoManager("Geometry", "default geometry");
163 }
164 if (fZ - Int_t(fZ) > 1E-3)
165 Warning("ctor", "Material %s defined with fractional Z=%f", GetName(), fZ);
166 if (GetElement())
167 GetElement()->SetUsed();
169}
170
171////////////////////////////////////////////////////////////////////////////////
172/// Constructor with state, temperature and pressure.
173///
174/// \param name material name.
175/// \param a atomic mass.
176/// \param z atomic number.
177/// \param rho material density in g/cm3.
178/// \param state
179/// \param temperature
180/// \param pressure
181
184 : TNamed(name, ""),
185 TAttFill(),
186 fIndex(0),
187 fA(a),
188 fZ(z),
189 fDensity(rho),
190 fRadLen(0.),
191 fIntLen(0.),
192 fTemperature(temperature),
193 fPressure(pressure),
194 fState(state),
195 fShader(nullptr),
196 fCerenkov(nullptr),
197 fElement(nullptr),
198 fUserExtension(nullptr),
199 fFWExtension(nullptr)
200{
201 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
202 fName = fName.Strip();
204 fIndex = -1;
205 SetRadLen(0, 0);
206 if (!gGeoManager) {
207 gGeoManager = new TGeoManager("Geometry", "default geometry");
208 }
209 if (fZ - Int_t(fZ) > 1E-3)
210 Warning("ctor", "Material %s defined with fractional Z=%f", GetName(), fZ);
211 if (GetElement())
212 GetElement()->SetUsed();
214}
215
216////////////////////////////////////////////////////////////////////////////////
217/// Constructor.
218///
219/// \param name material name.
220/// \param elem
221/// \param rho material density in g/cm3.
222
224 : TNamed(name, ""),
225 TAttFill(),
226 fIndex(0),
227 fA(0.),
228 fZ(0.),
229 fDensity(rho),
230 fRadLen(0.),
231 fIntLen(0.),
232 fTemperature(0.),
233 fPressure(0.),
234 fState(kMatStateUndefined),
235 fShader(nullptr),
236 fCerenkov(nullptr),
237 fElement(elem),
238 fUserExtension(nullptr),
239 fFWExtension(nullptr)
240{
241 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
242 fName = fName.Strip();
244 fIndex = -1;
245 fA = elem->A();
246 fZ = elem->Z();
247 SetRadLen(0, 0);
251 if (!gGeoManager) {
252 gGeoManager = new TGeoManager("Geometry", "default geometry");
253 }
254 if (fZ - Int_t(fZ) > 1E-3)
255 Warning("ctor", "Material %s defined with fractional Z=%f", GetName(), fZ);
256 if (GetElement())
257 GetElement()->SetUsed();
259}
260
261////////////////////////////////////////////////////////////////////////////////
262
264 : TNamed(gm),
265 TAttFill(gm),
266 fIndex(gm.fIndex),
267 fA(gm.fA),
268 fZ(gm.fZ),
269 fDensity(gm.fDensity),
270 fRadLen(gm.fRadLen),
271 fIntLen(gm.fIntLen),
272 fTemperature(gm.fTemperature),
273 fPressure(gm.fPressure),
274 fState(gm.fState),
275 fShader(gm.fShader),
276 fCerenkov(gm.fCerenkov),
277 fElement(gm.fElement),
278 fUserExtension(gm.fUserExtension->Grab()),
279 fFWExtension(gm.fFWExtension->Grab())
280
281{
282 // copy constructor
283 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
285 TIter next(&fProperties);
287 while ((property = (TNamed *)next()))
289}
290
291////////////////////////////////////////////////////////////////////////////////
292/// assignment operator
293
295{
296 if (this != &gm) {
298 TAttFill::operator=(gm);
299 fIndex = gm.fIndex;
300 fA = gm.fA;
301 fZ = gm.fZ;
302 fDensity = gm.fDensity;
303 fRadLen = gm.fRadLen;
304 fIntLen = gm.fIntLen;
305 fTemperature = gm.fTemperature;
306 fPressure = gm.fPressure;
307 fState = gm.fState;
308 fShader = gm.fShader;
309 fCerenkov = gm.fCerenkov;
310 fElement = gm.fElement;
311 fUserExtension = gm.fUserExtension->Grab();
312 fFWExtension = gm.fFWExtension->Grab();
314 TIter next(&fProperties);
316 while ((property = (TNamed *)next()))
318 }
319 return *this;
320}
321
322////////////////////////////////////////////////////////////////////////////////
323/// Destructor
324
326{
327 if (fUserExtension) {
329 fUserExtension = nullptr;
330 }
331 if (fFWExtension) {
333 fFWExtension = nullptr;
334 }
335}
336
337////////////////////////////////////////////////////////////////////////////////
338/// Connect user-defined extension to the material. The material "grabs" a copy, so
339/// the original object can be released by the producer. Release the previously
340/// connected extension if any.
341///
342/// NOTE: This interface is intended for user extensions and is guaranteed not
343/// to be used by TGeo
344
346{
347 if (fUserExtension)
349 fUserExtension = nullptr;
350 if (ext)
351 fUserExtension = ext->Grab();
352}
353
354//_____________________________________________________________________________
355const char *TGeoMaterial::GetPropertyRef(const char *property) const
356{
357 // Find reference for a given property
359 return (prop) ? prop->GetTitle() : nullptr;
360}
361
362//_____________________________________________________________________________
364{
365 // Find reference for a given property
367 if (!prop)
368 return nullptr;
369 return gGeoManager->GetGDMLMatrix(prop->GetTitle());
370}
371
372//_____________________________________________________________________________
374{
375 // Find reference for a given property
377 if (!prop)
378 return nullptr;
379 return gGeoManager->GetGDMLMatrix(prop->GetTitle());
380}
381
382//_____________________________________________________________________________
383const char *TGeoMaterial::GetConstPropertyRef(const char *property) const
384{
385 // Find reference for a given constant property
387 return (prop) ? prop->GetTitle() : nullptr;
388}
389
390//_____________________________________________________________________________
392{
393 // Find reference for a given constant property
395 if (!prop) {
396 if (err)
397 *err = kTRUE;
398 return 0.;
399 }
400 return gGeoManager->GetProperty(prop->GetTitle(), err);
401}
402
403//_____________________________________________________________________________
405{
406 // Find reference for a given constant property
408 if (!prop) {
409 if (err)
410 *err = kTRUE;
411 return 0.;
412 }
413 return gGeoManager->GetProperty(prop->GetTitle(), err);
414}
415
416//_____________________________________________________________________________
417bool TGeoMaterial::AddProperty(const char *property, const char *ref)
418{
421 Error("AddProperty", "Property %s already added to material %s", property, GetName());
422 return false;
423 }
425 return true;
426}
427
428//_____________________________________________________________________________
429bool TGeoMaterial::AddConstProperty(const char *property, const char *ref)
430{
433 Error("AddConstProperty", "Constant property %s already added to material %s", property, GetName());
434 return false;
435 }
437 return true;
438}
439
440////////////////////////////////////////////////////////////////////////////////
441/// Connect framework defined extension to the material. The material "grabs" a copy,
442/// so the original object can be released by the producer. Release the previously
443/// connected extension if any.
444///
445/// NOTE: This interface is intended for the use by TGeo and the users should
446/// NOT connect extensions using this method
447
449{
450 if (fFWExtension)
452 fFWExtension = nullptr;
453 if (ext)
454 fFWExtension = ext->Grab();
455}
456
457////////////////////////////////////////////////////////////////////////////////
458/// Get a copy of the user extension pointer. The user must call Release() on
459/// the copy pointer once this pointer is not needed anymore (equivalent to
460/// delete() after calling new())
461
463{
464 if (fUserExtension)
465 return fUserExtension->Grab();
466 return nullptr;
467}
468
469////////////////////////////////////////////////////////////////////////////////
470/// Get a copy of the framework extension pointer. The user must call Release() on
471/// the copy pointer once this pointer is not needed anymore (equivalent to
472/// delete() after calling new())
473
475{
476 if (fFWExtension)
477 return fFWExtension->Grab();
478 return nullptr;
479}
480
481////////////////////////////////////////////////////////////////////////////////
482/// Provide a pointer name containing uid.
483
485{
486 static TString name;
487 name.Form("pMat%d", GetUniqueID());
488 return name.Data();
489}
490
491////////////////////////////////////////////////////////////////////////////////
492/// Set radiation/absorption lengths. If the values are negative, their absolute value
493/// is taken, otherwise radlen is recomputed using G3 formula.
494
496{
499 // Check for vacuum
500 if (fA < 0.9 || fZ < 0.9) {
501 if (radlen < -1e5 || intlen < -1e-5) {
502 Error("SetRadLen", "Material %s: user values taken for vacuum: radlen=%g or intlen=%g - too small", GetName(),
504 return;
505 }
506 // Ignore positive values and take big numbers
507 if (radlen >= 0)
508 fRadLen = 1.E30;
509 if (intlen >= 0)
510 fIntLen = 1.E30;
511 return;
512 }
514 // compute radlen systematically with G3 formula for a valid material
515 if (radlen >= 0) {
516 // taken grom Geant3 routine GSMATE
517 constexpr Double_t alr2av = 1.39621E-03;
518 constexpr Double_t al183 = 5.20948;
521 // fRadLen is in TGeo units. Apply conversion factor in requested length-units
523 }
524 // Compute interaction length using the same formula as in GEANT4
525 if (intlen >= 0) {
526 constexpr Double_t lambda0 = 35. * TGeoUnit::g / TGeoUnit::cm2; // [g/cm^2]
527 Double_t nilinv = 0.0;
529 if (!elem) {
530 Fatal("SetRadLen", "Element not found for material %s", GetName());
531 return;
532 }
534 nilinv += nbAtomsPerVolume * TMath::Power(elem->Neff(), 0.6666667);
536 fIntLen = (nilinv <= 0) ? TGeoShape::Big() : (1.0 / nilinv);
537 // fIntLen is in TGeo units. Apply conversion factor in requested length-units
539 }
540}
541
542////////////////////////////////////////////////////////////////////////////////
543/// static function
544/// Compute Coulomb correction for pair production and Brem
545/// REFERENCE : EGS MANUAL SLAC 210 - UC32 - JUNE 78
546/// FORMULA 2.7.17
547
549{
552 Double_t az2 = az * az;
553 Double_t az4 = az2 * az2;
554 Double_t fp = (0.0083 * az4 + 0.20206 + 1. / (1. + az2)) * az2;
555 Double_t fm = (0.0020 * az4 + 0.0369) * az4;
556 return fp - fm;
557}
558
559////////////////////////////////////////////////////////////////////////////////
560/// return true if the other material has the same physical properties
561
563{
564 if (other == this)
565 return kTRUE;
566 if (other->IsMixture())
567 return kFALSE;
568 if (TMath::Abs(fA - other->GetA()) > 1E-3)
569 return kFALSE;
570 if (TMath::Abs(fZ - other->GetZ()) > 1E-3)
571 return kFALSE;
572 if (TMath::Abs(fDensity - other->GetDensity()) > 1E-6)
573 return kFALSE;
574 if (GetCerenkovProperties() != other->GetCerenkovProperties())
575 return kFALSE;
576 // if (fRadLen != other->GetRadLen()) return kFALSE;
577 // if (fIntLen != other->GetIntLen()) return kFALSE;
578 return kTRUE;
579}
580
581////////////////////////////////////////////////////////////////////////////////
582/// print characteristics of this material
583
584void TGeoMaterial::Print(const Option_t * /*option*/) const
585{
586 printf("Material %s %s A=%g Z=%g rho=%g radlen=%g intlen=%g index=%i\n", GetName(), GetTitle(), fA, fZ, fDensity,
588}
589
590////////////////////////////////////////////////////////////////////////////////
591/// Save a primitive as a C++ statement(s) on output stream "out".
592
593void TGeoMaterial::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
594{
596 return;
597 const char *name = GetPointerName();
598 out << "// Material: " << GetName() << std::endl;
599 out << " a = " << fA << ";" << std::endl;
600 out << " z = " << fZ << ";" << std::endl;
601 out << " density = " << fDensity << ";" << std::endl;
602 out << " radl = " << fRadLen << ";" << std::endl;
603 out << " absl = " << fIntLen << ";" << std::endl;
604
605 out << " auto " << name << " = new TGeoMaterial(\"" << GetName() << "\", a, z, density, radl, absl);" << std::endl;
606 out << " " << name << "->SetIndex(" << GetIndex() << ");" << std::endl;
608}
609
610////////////////////////////////////////////////////////////////////////////////
611/// Get some default color related to this material.
612
614{
615 Int_t id = 1 + gGeoManager->GetListOfMaterials()->IndexOf(this);
616 return (2 + id % 6);
617}
618
619////////////////////////////////////////////////////////////////////////////////
620/// Get a pointer to the element this material is made of.
621/// This second call is to avoid warnings to not call a virtual
622/// method from the constructor
623
625{
626 if (fElement)
627 return fElement;
629 return table->GetElement(Int_t(fZ));
630}
631
632////////////////////////////////////////////////////////////////////////////////
633/// Get a pointer to the element this material is made of.
634
636{
637 if (fElement)
638 return fElement;
640 return table->GetElement(Int_t(fZ));
641}
642
643////////////////////////////////////////////////////////////////////////////////
644/// Single interface to get element properties.
645
647{
648 a = fA;
649 z = fZ;
650 w = 1.;
651}
652
653////////////////////////////////////////////////////////////////////////////////
654/// Retrieve material index in the list of materials
655
657{
658 if (fIndex >= 0)
659 return fIndex;
661 fIndex = matlist->IndexOf(this);
662 return fIndex;
663}
664
665////////////////////////////////////////////////////////////////////////////////
666/// Create the material representing the decay product of this material at a
667/// given time. The precision represent the minimum cumulative branching ratio for
668/// which decay products are still taken into account.
669
671{
672 TObjArray *pop = new TObjArray();
673 if (!fElement || !fElement->IsRadioNuclide())
674 return this;
675 FillMaterialEvolution(pop, precision);
676 Int_t ncomp = pop->GetEntriesFast();
677 if (!ncomp)
678 return this;
680 Double_t *weight = new Double_t[ncomp];
681 Double_t amed = 0.;
682 Int_t i;
683 for (i = 0; i < ncomp; i++) {
684 el = (TGeoElementRN *)pop->At(i);
685 weight[i] = el->Ratio()->Concentration(time) * el->A();
686 amed += weight[i];
687 }
688 Double_t rho = fDensity * amed / fA;
689 TGeoMixture *mix = nullptr;
691 for (i = 0; i < ncomp; i++) {
692 if ((weight[i] / amed) < precision) {
693 amed -= weight[i];
694 ncomp1--;
695 }
696 }
697 if (ncomp1 < 2) {
698 el = (TGeoElementRN *)pop->At(0);
699 delete[] weight;
700 delete pop;
701 if (ncomp1 == 1)
702 return new TGeoMaterial(TString::Format("%s-evol", GetName()), el, rho);
703 return nullptr;
704 }
705 mix = new TGeoMixture(TString::Format("%s-evol", GetName()), ncomp, rho);
706 for (i = 0; i < ncomp; i++) {
707 weight[i] /= amed;
708 if (weight[i] < precision)
709 continue;
710 el = (TGeoElementRN *)pop->At(i);
711 mix->AddElement(el, weight[i]);
712 }
713 delete[] weight;
714 delete pop;
715 return mix;
716}
717
718////////////////////////////////////////////////////////////////////////////////
719/// Fills a user array with all the elements deriving from the possible
720/// decay of the top element composing the mixture. Each element contained
721/// by `<population>` may be a radionuclide having a Bateman solution attached.
722/// The precision represent the minimum cumulative branching ratio for
723/// which decay products are still taken into account.
724/// To visualize the time evolution of each decay product one can use:
725/// ~~~ {.cpp}
726/// TGeoElement *elem = population->At(index);
727/// TGeoElementRN *elemrn = 0;
728/// if (elem->IsRadioNuclide()) elemrn = (TGeoElementRN*)elem;
729/// ~~~
730/// One can get Ni/N1(t=0) at any moment of time. Ni is the number of atoms
731/// of one of the decay products, N1(0) is the number of atoms of the top
732/// element at t=0.
733/// ~~~ {.cpp}
734/// Double_t fraction_weight = elemrn->Ratio()->Concentration(time);
735/// ~~~
736/// One can also display the time evolution of the fractional weight:
737/// ~~~ {.cpp}
738/// elemrn->Ratio()->Draw(option);
739/// ~~~
740
742{
743 if (population->GetEntriesFast()) {
744 Error("FillMaterialEvolution", "Provide an empty array !");
745 return;
746 }
750 TIter next(table->GetElementsRN());
751 while ((elemrn = (TGeoElementRN *)next()))
752 elemrn->ResetRatio();
753 elem = GetElement();
754 if (!elem) {
755 Fatal("FillMaterialEvolution", "Element not found for material %s", GetName());
756 return;
757 }
758 if (!elem->IsRadioNuclide()) {
759 population->Add(elem);
760 return;
761 }
763 elemrn->FillPopulation(population, precision);
764}
765
766/** \class TGeoMixture
767\ingroup Materials_classes
768
769Mixtures of elements.
770
771*/
772
773
774////////////////////////////////////////////////////////////////////////////////
775/// Default constructor
776
778{
779 fNelements = 0;
780 fZmixture = nullptr;
781 fAmixture = nullptr;
782 fWeights = nullptr;
783 fNatoms = nullptr;
784 fVecNbOfAtomsPerVolume = nullptr;
785 fElements = nullptr;
786}
787
788////////////////////////////////////////////////////////////////////////////////
789/// constructor
790
792{
793 fZmixture = nullptr;
794 fAmixture = nullptr;
795 fWeights = nullptr;
796 fNelements = 0;
797 fNatoms = nullptr;
798 fVecNbOfAtomsPerVolume = nullptr;
799 fDensity = rho;
800 fElements = nullptr;
801 if (fDensity < 0)
802 fDensity = 0.001;
803}
804
805////////////////////////////////////////////////////////////////////////////////
806/// Destructor
807
809{
810 if (fZmixture)
811 delete[] fZmixture;
812 if (fAmixture)
813 delete[] fAmixture;
814 if (fWeights)
815 delete[] fWeights;
816 if (fNatoms)
817 delete[] fNatoms;
819 delete[] fVecNbOfAtomsPerVolume;
820 if (fElements)
821 delete fElements;
822}
823
824////////////////////////////////////////////////////////////////////////////////
825/// Compute effective A/Z and radiation length
826
828{
829 constexpr const Double_t na = TGeoUnit::Avogadro;
830 constexpr const Double_t alr2av = 1.39621E-03;
831 constexpr const Double_t al183 = 5.20948;
832 constexpr const Double_t lambda0 = 35. * TGeoUnit::g / TGeoUnit::cm2; // [g/cm^2]
833 Double_t radinv = 0.0;
834 Double_t nilinv = 0.0;
836 fA = 0;
837 fZ = 0;
838 for (Int_t j = 0; j < fNelements; j++) {
839 if (fWeights[j] <= 0)
840 continue;
841 fA += fWeights[j] * fAmixture[j];
842 fZ += fWeights[j] * fZmixture[j];
844 nilinv += nbAtomsPerVolume * TMath::Power(GetElement(j)->Neff(), 0.6666667);
846 Double_t alz = TMath::Log(zc) / 3.;
847 Double_t xinv =
849 radinv += xinv * fWeights[j];
850 }
852 fRadLen = (radinv <= 0) ? TGeoShape::Big() : 1.0 / radinv;
853 // fRadLen is in TGeo units. Apply conversion factor in requested length-units
855
856 // Compute interaction length
858 fIntLen = (nilinv <= 0) ? TGeoShape::Big() : 1.0 / nilinv;
859 // fIntLen is in TGeo units. Apply conversion factor in requested length-units
861}
862
863////////////////////////////////////////////////////////////////////////////////
864/// add an element to the mixture using fraction by weight
865/// Check if the element is already defined
866
868{
870
871 // Check preconditions
872 if (weight < 0e0) {
873 Fatal("AddElement", "Cannot add element with negative weight %g to mixture %s", weight, GetName());
874 } else if (weight < std::numeric_limits<Double_t>::epsilon()) {
875 return;
876 } else if (z < 1 || z > table->GetNelements() - 1) {
877 Fatal("AddElement", "Cannot add element having Z=%d to mixture %s", (Int_t)z, GetName());
878 }
879 Int_t i;
880 for (i = 0; i < fNelements; i++) {
881 if (!fElements && TMath::Abs(z - fZmixture[i]) < 1.e-6 && TMath::Abs(a - fAmixture[i]) < 1.e-6) {
882 fWeights[i] += weight;
884 return;
885 }
886 }
887 if (!fNelements) {
888 fZmixture = new Double_t[1];
889 fAmixture = new Double_t[1];
890 fWeights = new Double_t[1];
891 } else {
895 Double_t *weights = new Double_t[nelements];
896 for (Int_t j = 0; j < fNelements; j++) {
897 zmixture[j] = fZmixture[j];
898 amixture[j] = fAmixture[j];
899 weights[j] = fWeights[j];
900 }
901 delete[] fZmixture;
902 delete[] fAmixture;
903 delete[] fWeights;
906 fWeights = weights;
907 }
908
909 fNelements++;
910 i = fNelements - 1;
911 fZmixture[i] = z;
912 fAmixture[i] = a;
913 fWeights[i] = weight;
914 if (z - Int_t(z) > 1E-3)
915 Warning("DefineElement", "Mixture %s has element defined with fractional Z=%f", GetName(), z);
917 table->GetElement((Int_t)z)->SetDefined();
918
919 // compute equivalent radiation length (taken from Geant3/GSMIXT)
921}
922
923////////////////////////////////////////////////////////////////////////////////
924/// Define one component of the mixture as an existing material/mixture.
925
927{
929 Double_t a, z;
930
931 // Check preconditions
932 if (!mat) {
933 Fatal("AddElement", "Cannot add INVALID material to mixture %s", GetName());
934 } else if (weight < 0e0) {
935 Fatal("AddElement", "Cannot add material %s with negative weight %g to mixture %s", mat->GetName(), weight,
936 GetName());
937 } else if (weight < std::numeric_limits<Double_t>::epsilon()) {
938 return;
939 }
940 if (!mat->IsMixture()) {
941 elem = mat->GetBaseElement();
942 if (elem) {
943 AddElement(elem, weight);
944 } else {
945 a = mat->GetA();
946 z = mat->GetZ();
947 AddElement(a, z, weight);
948 }
949 return;
950 }
951 // The material is a mixture.
952 TGeoMixture *mix = (TGeoMixture *)mat;
954 Int_t nelem = mix->GetNelements();
956 Int_t i, j;
957 // loop the elements of the daughter mixture
958 for (i = 0; i < nelem; i++) {
959 elfound = kFALSE;
960 elnew = mix->GetElement(i);
961 if (!elnew)
962 continue;
963 // check if we have the element already defined in the parent mixture
964 for (j = 0; j < fNelements; j++) {
965 if (fWeights[j] < 0e0)
966 continue;
967 elem = GetElement(j);
968 if (elem == elnew) {
969 // element found, compute new weight
970 fWeights[j] += weight * (mix->GetWmixt())[i];
971 elfound = kTRUE;
972 break;
973 }
974 }
975 if (elfound)
976 continue;
977 // element not found, define it
978 wnew = weight * (mix->GetWmixt())[i];
980 }
981}
982
983////////////////////////////////////////////////////////////////////////////////
984/// add an element to the mixture using fraction by weight
985
987{
990 if (!fElements)
991 fElements = new TObjArray(128);
993
994 // Check preconditions
995 if (!elem) {
996 Fatal("AddElement", "Cannot add INVALID element to mixture %s", GetName());
997 } else if (weight < 0e0) {
998 Fatal("AddElement", "Cannot add element %s with negative weight %g to mixture %s", elem->GetName(), weight,
999 GetName());
1000 } else if (weight < std::numeric_limits<Double_t>::epsilon()) {
1001 return;
1002 }
1003 // If previous elements were defined by A/Z, add corresponding TGeoElements
1004 for (Int_t i = 0; i < fNelements; i++) {
1006 if (!elemold) {
1007 // Add element with corresponding Z in the list
1008 fElements->AddAt(elemold = table->GetElement((Int_t)fZmixture[i]), i);
1009 elemold->SetDefined();
1010 }
1011 if (elemold == elem) {
1012 fWeights[i] += weight;
1013 exist = kTRUE;
1014 }
1015 }
1016 if (!exist) {
1018 AddElement(elem->A(), elem->Z(), weight);
1019 } else {
1021 }
1022}
1023
1024////////////////////////////////////////////////////////////////////////////////
1025/// Add a mixture element by number of atoms in the chemical formula.
1026
1028{
1029 Int_t i, j;
1030 Double_t amol;
1033 if (!fElements)
1034 fElements = new TObjArray(128);
1035 // Check if the element is already defined
1036 for (i = 0; i < fNelements; i++) {
1038 if (!elemold)
1039 fElements->AddAt(table->GetElement((Int_t)fZmixture[i]), i);
1040 else if (elemold != elem)
1041 continue;
1042 if ((elem == elemold) ||
1043 (TMath::Abs(elem->Z() - fZmixture[i]) < 1.e-6 && TMath::Abs(elem->A() - fAmixture[i]) < 1.e-6)) {
1044 fNatoms[i] += natoms;
1045 amol = 0.;
1046 for (j = 0; j < fNelements; j++)
1047 amol += fAmixture[j] * fNatoms[j];
1048 for (j = 0; j < fNelements; j++)
1049 fWeights[j] = fNatoms[j] * fAmixture[j] / amol;
1051 return;
1052 }
1053 }
1054 // New element
1055 if (!fNelements) {
1056 fZmixture = new Double_t[1];
1057 fAmixture = new Double_t[1];
1058 fWeights = new Double_t[1];
1059 fNatoms = new Int_t[1];
1060 } else {
1061 if (!fNatoms) {
1062 Fatal("AddElement", "Cannot add element by natoms in mixture %s after defining elements by weight", GetName());
1063 return;
1064 }
1068 Double_t *weights = new Double_t[nelements];
1069 Int_t *nnatoms = new Int_t[nelements];
1070 for (j = 0; j < fNelements; j++) {
1071 zmixture[j] = fZmixture[j];
1072 amixture[j] = fAmixture[j];
1073 weights[j] = fWeights[j];
1074 nnatoms[j] = fNatoms[j];
1075 }
1076 delete[] fZmixture;
1077 delete[] fAmixture;
1078 delete[] fWeights;
1079 delete[] fNatoms;
1082 fWeights = weights;
1083 fNatoms = nnatoms;
1084 }
1085 fNelements++;
1086 Int_t iel = fNelements - 1;
1087 fZmixture[iel] = elem->Z();
1088 fAmixture[iel] = elem->A();
1089 fNatoms[iel] = natoms;
1091 amol = 0.;
1092 for (i = 0; i < fNelements; i++) {
1093 if (fNatoms[i] <= 0)
1094 return;
1095 amol += fAmixture[i] * fNatoms[i];
1096 }
1097 for (i = 0; i < fNelements; i++)
1098 fWeights[i] = fNatoms[i] * fAmixture[i] / amol;
1099 table->GetElement(elem->Z())->SetDefined();
1101}
1102
1103////////////////////////////////////////////////////////////////////////////////
1104/// Define the mixture element at index iel by number of atoms in the chemical formula.
1105
1107{
1109 TGeoElement *elem = table->GetElement(z);
1110 if (!elem) {
1111 Fatal("DefineElement", "In mixture %s, element with Z=%i not found", GetName(), z);
1112 return;
1113 }
1115}
1116
1117////////////////////////////////////////////////////////////////////////////////
1118/// Retrieve the pointer to the element corresponding to component I.
1119
1121{
1122 if (i < 0 || i >= fNelements) {
1123 Error("GetElement", "Mixture %s has only %d elements", GetName(), fNelements);
1124 return nullptr;
1125 }
1126 TGeoElement *elem = nullptr;
1127 if (fElements)
1128 elem = (TGeoElement *)fElements->At(i);
1129 if (elem)
1130 return elem;
1132 return table->GetElement(Int_t(fZmixture[i]));
1133}
1134
1135////////////////////////////////////////////////////////////////////////////////
1136/// Get specific activity (in Bq/gram) for the whole mixture (no argument) or
1137/// for a given component.
1138
1140{
1141 if (i >= 0 && i < fNelements)
1142 return fWeights[i] * GetElement(i)->GetSpecificActivity();
1143 Double_t sa = 0;
1144 for (Int_t iel = 0; iel < fNelements; iel++) {
1146 }
1147 return sa;
1148}
1149
1150////////////////////////////////////////////////////////////////////////////////
1151/// Return true if the other material has the same physical properties
1152
1154{
1155 if (other->IsEqual(this))
1156 return kTRUE;
1157 if (!other->IsMixture())
1158 return kFALSE;
1159 TGeoMixture *mix = (TGeoMixture *)other;
1160 if (!mix)
1161 return kFALSE;
1162 if (fNelements != mix->GetNelements())
1163 return kFALSE;
1164 if (TMath::Abs(fA - other->GetA()) > 1E-3)
1165 return kFALSE;
1166 if (TMath::Abs(fZ - other->GetZ()) > 1E-3)
1167 return kFALSE;
1168 if (TMath::Abs(fDensity - other->GetDensity()) > 1E-6)
1169 return kFALSE;
1170 if (GetCerenkovProperties() != other->GetCerenkovProperties())
1171 return kFALSE;
1172 // if (fRadLen != other->GetRadLen()) return kFALSE;
1173 // if (fIntLen != other->GetIntLen()) return kFALSE;
1174 for (Int_t i = 0; i < fNelements; i++) {
1175 if (TMath::Abs(fZmixture[i] - (mix->GetZmixt())[i]) > 1E-3)
1176 return kFALSE;
1177 if (TMath::Abs(fAmixture[i] - (mix->GetAmixt())[i]) > 1E-3)
1178 return kFALSE;
1179 if (TMath::Abs(fWeights[i] - (mix->GetWmixt())[i]) > 1E-3)
1180 return kFALSE;
1181 }
1182 return kTRUE;
1183}
1184
1185////////////////////////////////////////////////////////////////////////////////
1186/// print characteristics of this material
1187
1188void TGeoMixture::Print(const Option_t * /*option*/) const
1189{
1190 printf("Mixture %s %s Aeff=%g Zeff=%g rho=%g radlen=%g intlen=%g index=%i\n", GetName(), GetTitle(), fA, fZ,
1192 for (Int_t i = 0; i < fNelements; i++) {
1193 if (fElements && fElements->At(i)) {
1194 printf(" Element #%i : %s Z=%6.2f A=%6.2f w=%6.3f\n", i, GetElement(i)->GetName(), fZmixture[i],
1195 fAmixture[i], fWeights[i]);
1196 continue;
1197 }
1198 if (fNatoms)
1199 printf(" Element #%i : %s Z=%6.2f A=%6.2f w=%6.3f natoms=%d\n", i, GetElement(i)->GetName(), fZmixture[i],
1200 fAmixture[i], fWeights[i], fNatoms[i]);
1201 else
1202 printf(" Element #%i : %s Z=%6.2f A=%6.2f w=%6.3f\n", i, GetElement(i)->GetName(), fZmixture[i],
1203 fAmixture[i], fWeights[i]);
1204 }
1205}
1206
1207////////////////////////////////////////////////////////////////////////////////
1208/// Save a primitive as a C++ statement(s) on output stream "out".
1209
1210void TGeoMixture::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
1211{
1213 return;
1214 const char *name = GetPointerName();
1215 out << "// Mixture: " << GetName() << std::endl;
1216 out << " nel = " << fNelements << ";" << std::endl;
1217 out << " density = " << fDensity << ";" << std::endl;
1218 out << " auto " << name << " = new TGeoMixture(\"" << GetName() << "\", nel, density);" << std::endl;
1219 for (Int_t i = 0; i < fNelements; i++) {
1221 out << " a = " << fAmixture[i] << "; z = " << fZmixture[i] << "; w = " << fWeights[i] << "; // "
1222 << el->GetName() << std::endl;
1223 out << " " << name << "->DefineElement(" << i << ",a,z,w);" << std::endl;
1224 }
1225 out << " " << name << "->SetIndex(" << GetIndex() << ");" << std::endl;
1227}
1228
1229////////////////////////////////////////////////////////////////////////////////
1230/// Create the mixture representing the decay product of this material at a
1231/// given time. The precision represent the minimum cumulative branching ratio for
1232/// which decay products are still taken into account.
1233
1235{
1236 TObjArray *pop = new TObjArray();
1237 FillMaterialEvolution(pop, precision);
1238 Int_t ncomp = pop->GetEntriesFast();
1239 if (!ncomp)
1240 return this;
1243 Double_t *weight = new Double_t[ncomp];
1244 Double_t amed = 0.;
1245 Int_t i, j;
1246 for (i = 0; i < ncomp; i++) {
1247 elem = (TGeoElement *)pop->At(i);
1248 if (!elem->IsRadioNuclide()) {
1249 j = fElements->IndexOf(elem);
1250 weight[i] = fWeights[j] * fAmixture[0] / fWeights[0];
1251 } else {
1252 el = (TGeoElementRN *)elem;
1253 weight[i] = el->Ratio()->Concentration(time) * el->A();
1254 }
1255 amed += weight[i];
1256 }
1257 Double_t rho = fDensity * fWeights[0] * amed / fAmixture[0];
1258 TGeoMixture *mix = nullptr;
1259 Int_t ncomp1 = ncomp;
1260 for (i = 0; i < ncomp; i++) {
1261 if ((weight[i] / amed) < precision) {
1262 amed -= weight[i];
1263 ncomp1--;
1264 }
1265 }
1266 if (ncomp1 < 2) {
1267 el = (TGeoElementRN *)pop->At(0);
1268 delete[] weight;
1269 delete pop;
1270 if (ncomp1 == 1)
1271 return new TGeoMaterial(TString::Format("%s-evol", GetName()), el, rho);
1272 return nullptr;
1273 }
1274 mix = new TGeoMixture(TString::Format("%s-evol", GetName()), ncomp, rho);
1275 for (i = 0; i < ncomp; i++) {
1276 weight[i] /= amed;
1277 if (weight[i] < precision)
1278 continue;
1279 el = (TGeoElementRN *)pop->At(i);
1280 mix->AddElement(el, weight[i]);
1281 }
1282 delete[] weight;
1283 delete pop;
1284 return mix;
1285}
1286
1287////////////////////////////////////////////////////////////////////////////////
1288/// Fills a user array with all the elements deriving from the possible
1289/// decay of the top elements composing the mixture. Each element contained
1290/// by `<population>` may be a radionuclide having a Bateman solution attached.
1291/// The precision represent the minimum cumulative branching ratio for
1292/// which decay products are still taken into account.
1293/// To visualize the time evolution of each decay product one can use:
1294/// ~~~ {.cpp}
1295/// TGeoElement *elem = population->At(index);
1296/// TGeoElementRN *elemrn = 0;
1297/// if (elem->IsRadioNuclide()) elemrn = (TGeoElementRN*)elem;
1298/// ~~~
1299/// One can get Ni/N1(t=0) at any moment of time. Ni is the number of atoms
1300/// of one of the decay products, N1(0) is the number of atoms of the first top
1301/// element at t=0.
1302/// ~~~ {.cpp}
1303/// Double_t fraction_weight = elemrn->Ratio()->Concentration(time);
1304/// ~~~
1305/// One can also display the time evolution of the fractional weight:
1306/// ~~~ {.cpp}
1307/// elemrn->Ratio()->Draw(option);
1308/// ~~~
1309
1311{
1312 if (population->GetEntriesFast()) {
1313 Error("FillMaterialEvolution", "Provide an empty array !");
1314 return;
1315 }
1319 TIter next(table->GetElementsRN());
1320 while ((elemrn = (TGeoElementRN *)next()))
1321 elemrn->ResetRatio();
1322 Double_t factor;
1323 for (Int_t i = 0; i < fNelements; i++) {
1324 elem = GetElement(i);
1325 if (!elem->IsRadioNuclide()) {
1326 population->Add(elem);
1327 continue;
1328 }
1330 factor = fWeights[i] * fAmixture[0] / (fWeights[0] * fAmixture[i]);
1331 elemrn->FillPopulation(population, precision, factor);
1332 }
1333}
1334
1335////////////////////////////////////////////////////////////////////////////////
1336/// static function
1337/// Compute screening factor for pair production and Bremsstrahlung
1338/// REFERENCE : EGS MANUAL SLAC 210 - UC32 - JUNE 78
1339/// FORMULA 2.7.22
1340
1342{
1343 const Double_t al183 = 5.20948, al1440 = 7.27239;
1344 Double_t alz = TMath::Log(z) / 3.;
1345 Double_t factor = (al1440 - 2 * alz) / (al183 - alz - TGeoMaterial::Coulomb(z));
1346 return factor;
1347}
1348
1349////////////////////////////////////////////////////////////////////////////////
1350/// Compute Derived Quantities as in Geant4
1351
1353{
1354 const Double_t Na =
1356
1358 delete[] fVecNbOfAtomsPerVolume;
1359
1361
1362 // Formula taken from G4Material.cxx L312
1363 double sumweights = 0;
1364 for (Int_t i = 0; i < fNelements; ++i) {
1365 sumweights += fWeights[i];
1366 fVecNbOfAtomsPerVolume[i] = Na * fDensity * fWeights[i] / ((TGeoElement *)fElements->At(i))->A();
1367 }
1368 if (TMath::Abs(sumweights - 1) > 0.001)
1369 Warning("ComputeDerivedQuantities", "Mixture %s: sum of weights is: %g", GetName(), sumweights);
1372}
1373
1374////////////////////////////////////////////////////////////////////////////////
1375/// Compute Radiation Length based on Geant4 formula
1376
1378{
1379 // Formula taken from G4Material.cxx L556
1380 Double_t radinv = 0.0;
1381 // GetfRadTsai is in units of cm2 due to <unit>::alpha_rcl2. Correction must be applied to end up in TGeo cm.
1383 for (Int_t i = 0; i < fNelements; ++i) {
1384 radinv += fVecNbOfAtomsPerVolume[i] * ((TGeoElement *)fElements->At(i))->GetfRadTsai() / denom;
1385 }
1386 fRadLen = (radinv <= 0.0 ? DBL_MAX : 1.0 / radinv);
1387 // fRadLen is in TGeo units. Apply conversion factor in requested length-units
1389}
1390
1391////////////////////////////////////////////////////////////////////////////////
1392/// Compute Nuclear Interaction Length based on Geant4 formula
1394{
1395 // Formula taken from G4Material.cxx L567
1396 constexpr Double_t lambda0 = 35. * TGeoUnit::g / TGeoUnit::cm2; // [g/cm^2]
1397 const Double_t twothird = 2.0 / 3.0;
1398 Double_t NILinv = 0.0;
1399 for (Int_t i = 0; i < fNelements; ++i) {
1400 Int_t Z = static_cast<Int_t>(((TGeoElement *)fElements->At(i))->Z() + 0.5);
1401 Double_t A = ((TGeoElement *)fElements->At(i))->Neff();
1402 if (1 == Z) {
1404 } else {
1406 }
1407 }
1409 fIntLen = (NILinv <= 0.0 ? DBL_MAX : 1.0 / NILinv);
1410 // fIntLen is in TGeo units. Apply conversion factor in requested length-units
1412}
#define a(i)
Definition RSha256.hxx:99
#define e(i)
Definition RSha256.hxx:103
int Int_t
Signed integer 4 bytes (int)
Definition RtypesCore.h:59
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
constexpr Bool_t kTRUE
Definition RtypesCore.h:107
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h prop
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t Atom_t Atom_t Time_t property
char name[80]
Definition TGX11.cxx:110
R__EXTERN TGeoManager * gGeoManager
static const Double_t STP_temperature
static const Double_t STP_pressure
Fill Area Attributes class.
Definition TAttFill.h:20
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
This class is used in the process of reading and writing the GDML "matrix" tag.
Definition TGDMLMatrix.h:33
a radionuclide
table of elements
TGeoElement * GetElement(Int_t z)
TObjArray * GetElementsRN() const
Int_t GetNelements() const
Base class for chemical elements.
Definition TGeoElement.h:31
virtual Double_t GetSpecificActivity() const
Definition TGeoElement.h:74
Double_t A() const
Definition TGeoElement.h:66
void SetDefined(Bool_t flag=kTRUE)
Definition TGeoElement.h:80
virtual Bool_t IsRadioNuclide() const
Definition TGeoElement.h:77
void SetUsed(Bool_t flag=kTRUE)
Definition TGeoElement.h:81
ABC for user objects attached to TGeoVolume or TGeoNode.
virtual TGeoExtension * Grab()=0
virtual void Release() const =0
The manager class for any TGeo geometry.
Definition TGeoManager.h:45
static EDefaultUnits GetDefaultUnits()
TGeoElementTable * GetElementTable()
Returns material table. Creates it if not existing.
static void SetDefaultUnits(EDefaultUnits new_value)
Int_t AddMaterial(const TGeoMaterial *material)
Add a material to the list. Returns index of the material in list.
Double_t GetProperty(const char *name, Bool_t *error=nullptr) const
Get a user-defined property.
TGDMLMatrix * GetGDMLMatrix(const char *name) const
Get GDML matrix with a given name;.
TList * GetListOfMaterials() const
Base class describing materials.
Double_t GetConstProperty(const char *property, Bool_t *error=nullptr) const
void SetUserExtension(TGeoExtension *ext)
Connect user-defined extension to the material.
virtual TObject * GetCerenkovProperties() const
EGeoMaterialState fState
const char * GetPointerName() const
Provide a pointer name containing uid.
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save a primitive as a C++ statement(s) on output stream "out".
static Double_t ScreenFactor(Double_t z)
static function Compute screening factor for pair production and Bremsstrahlung REFERENCE : EGS MANUA...
const char * GetConstPropertyRef(const char *property) const
Double_t fPressure
void SetFWExtension(TGeoExtension *ext)
Connect framework defined extension to the material.
bool AddConstProperty(const char *property, const char *ref)
Double_t fTemperature
virtual void GetElementProp(Double_t &a, Double_t &z, Double_t &w, Int_t i=0)
Single interface to get element properties.
virtual TGeoMaterial * DecayMaterial(Double_t time, Double_t precision=0.001)
Create the material representing the decay product of this material at a given time.
TList fProperties
bool AddProperty(const char *property, const char *ref)
Double_t fZ
void SetRadLen(Double_t radlen, Double_t intlen=0.)
Set radiation/absorption lengths.
TGeoElement * GetElement() const
Get a pointer to the element this material is made of.
virtual Bool_t IsEq(const TGeoMaterial *other) const
return true if the other material has the same physical properties
const char * GetPropertyRef(const char *property) const
Double_t fA
TGDMLMatrix * GetProperty(const char *name) const
TObject * fCerenkov
Double_t fDensity
void SetUsed(Bool_t flag=kTRUE)
virtual void FillMaterialEvolution(TObjArray *population, Double_t precision=0.001)
Fills a user array with all the elements deriving from the possible decay of the top element composin...
Int_t GetIndex()
Retrieve material index in the list of materials.
Double_t fIntLen
TGeoExtension * fUserExtension
TList fConstProperties
TGeoElement * fElement
TGeoMaterial & operator=(const TGeoMaterial &)
assignment operator
TObject * fShader
TGeoExtension * GrabUserExtension() const
Get a copy of the user extension pointer.
Double_t fRadLen
TGeoExtension * GrabFWExtension() const
Get a copy of the framework extension pointer.
virtual Int_t GetDefaultColor() const
Get some default color related to this material.
TGeoMaterial()
Default constructor.
static Double_t Coulomb(Double_t z)
static function Compute Coulomb correction for pair production and Brem REFERENCE : EGS MANUAL SLAC 2...
TGeoExtension * fFWExtension
Transient user-defined extension to materials.
void Print(const Option_t *option="") const override
print characteristics of this material
~TGeoMaterial() override
Destructor.
Mixtures of elements.
TObjArray * fElements
void ComputeNuclearInterLength()
Compute Nuclear Interaction Length based on Geant4 formula.
~TGeoMixture() override
Destructor.
Double_t * GetZmixt() const
Bool_t IsEq(const TGeoMaterial *other) const override
Return true if the other material has the same physical properties.
void AddElement(Double_t a, Double_t z, Double_t weight)
add an element to the mixture using fraction by weight Check if the element is already defined
Double_t * fZmixture
TGeoMixture()
Default constructor.
void ComputeDerivedQuantities()
Compute Derived Quantities as in Geant4.
Int_t GetNelements() const override
TGeoMaterial * DecayMaterial(Double_t time, Double_t precision=0.001) override
Create the mixture representing the decay product of this material at a given time.
void ComputeRadiationLength()
Compute Radiation Length based on Geant4 formula.
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save a primitive as a C++ statement(s) on output stream "out".
Double_t * fVecNbOfAtomsPerVolume
Double_t * fAmixture
void FillMaterialEvolution(TObjArray *population, Double_t precision=0.001) override
Fills a user array with all the elements deriving from the possible decay of the top elements composi...
void AverageProperties()
Compute effective A/Z and radiation length.
void Print(const Option_t *option="") const override
print characteristics of this material
Double_t * GetWmixt() const
Double_t GetSpecificActivity(Int_t i=-1) const override
Get specific activity (in Bq/gram) for the whole mixture (no argument) or for a given component.
Int_t * fNatoms
Double_t * fWeights
void DefineElement(Int_t iel, Double_t a, Double_t z, Double_t weight)
Double_t * GetAmixt() const
TGeoElement * GetElement(Int_t i=0) const override
Retrieve the pointer to the element corresponding to component I.
static Double_t Big()
Definition TGeoShape.h:94
A doubly linked list.
Definition TList.h:38
TObject * FindObject(const char *name) const override
Find an object in this list using its name.
Definition TList.cxx:575
void Add(TObject *obj) override
Definition TList.h:81
TObject * At(Int_t idx) const override
Returns the object at position idx. Returns 0 if idx is out of range.
Definition TList.cxx:354
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
const char * GetName() const override
Returns name of object.
Definition TNamed.h:49
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:50
TNamed()
Definition TNamed.h:38
TString fName
Definition TNamed.h:32
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition TNamed.cxx:50
An array of TObjects.
Definition TObjArray.h:31
Int_t GetEntriesFast() const
Definition TObjArray.h:58
Int_t IndexOf(const TObject *obj) const override
void AddAt(TObject *obj, Int_t idx) override
Add object at position ids.
virtual void AddAtAndExpand(TObject *obj, Int_t idx)
Add object at position idx.
TObject * At(Int_t idx) const override
Definition TObjArray.h:164
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition TObject.h:202
virtual UInt_t GetUniqueID() const
Return the unique object id.
Definition TObject.cxx:475
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:1057
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition TObject.cxx:864
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
virtual void Fatal(const char *method, const char *msgfmt,...) const
Issue fatal error message.
Definition TObject.cxx:1099
virtual Int_t IndexOf(const TObject *obj) const
Return index of object in collection.
Basic string class.
Definition TString.h:138
TSubString Strip(EStripType s=kTrailing, char c=' ') const
Return a substring of self stripped at beginning and/or end.
Definition TString.cxx:1170
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:2384
static constexpr double fine_structure_const
static constexpr double cm2
static constexpr double Avogadro
static constexpr double cm
static constexpr double amu
static constexpr double cm
static constexpr double fine_structure_const
static constexpr double cm2
static constexpr double Avogadro
static constexpr double g
Double_t Exp(Double_t x)
Returns the base-e exponential function of x, which is e raised to the power x.
Definition TMath.h:720
Double_t Log(Double_t x)
Returns the natural logarithm of x.
Definition TMath.h:767
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Returns x raised to the power y.
Definition TMath.h:732
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:124