Logo ROOT  
Reference Guide
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
56
57////////////////////////////////////////////////////////////////////////////////
58/// Default constructor.
59
61 :TNamed(), 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(NULL),
72 fCerenkov(NULL),
73 fElement(NULL),
74 fUserExtension(0),
75 fFWExtension(0)
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, ""), TAttFill(),
92 fIndex(0),
93 fA(0.),
94 fZ(0.),
95 fDensity(0.),
96 fRadLen(0.),
97 fIntLen(0.),
98 fTemperature(0.),
99 fPressure(0.),
100 fState(kMatStateUndefined),
101 fShader(NULL),
102 fCerenkov(NULL),
103 fElement(NULL),
104 fUserExtension(0),
105 fFWExtension(0)
106{
107 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
108 fName = fName.Strip();
110 fIndex = -1;
114
115 if (!gGeoManager) {
116 gGeoManager = new TGeoManager("Geometry", "default geometry");
117 }
119}
120
121////////////////////////////////////////////////////////////////////////////////
122/// Constructor.
123///
124/// \param name material name.
125/// \param a atomic mass.
126/// \param z atomic number.
127/// \param rho material density in g/cm3.
128/// \param radlen
129/// \param intlen
130
132 Double_t rho, Double_t radlen, Double_t intlen)
133 :TNamed(name, ""), TAttFill(),
134 fIndex(0),
135 fA(a),
136 fZ(z),
137 fDensity(rho),
138 fRadLen(0.),
139 fIntLen(0.),
140 fTemperature(0.),
141 fPressure(0.),
142 fState(kMatStateUndefined),
143 fShader(NULL),
144 fCerenkov(NULL),
145 fElement(NULL),
146 fUserExtension(0),
147 fFWExtension(0)
148{
149 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
150 fName = fName.Strip();
152 fIndex = -1;
153 fA = a;
154 fZ = z;
155 fDensity = rho;
159 SetRadLen(radlen, intlen);
160 if (!gGeoManager) {
161 gGeoManager = new TGeoManager("Geometry", "default geometry");
162 }
163 if (fZ - Int_t(fZ) > 1E-3)
164 Warning("ctor", "Material %s defined with fractional Z=%f", GetName(), fZ);
165 if (GetElement()) GetElement()->SetUsed();
167}
168
169////////////////////////////////////////////////////////////////////////////////
170/// Constructor with state, temperature and pressure.
171///
172/// \param name material name.
173/// \param a atomic mass.
174/// \param z atomic number.
175/// \param rho material density in g/cm3.
176/// \param state
177/// \param temperature
178/// \param pressure
179
180
182 EGeoMaterialState state, Double_t temperature, Double_t pressure)
183 :TNamed(name, ""), TAttFill(),
184 fIndex(0),
185 fA(a),
186 fZ(z),
187 fDensity(rho),
188 fRadLen(0.),
189 fIntLen(0.),
190 fTemperature(temperature),
191 fPressure(pressure),
192 fState(state),
193 fShader(NULL),
194 fCerenkov(NULL),
195 fElement(NULL),
196 fUserExtension(0),
197 fFWExtension(0)
198{
199 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
200 fName = fName.Strip();
202 fIndex = -1;
203 SetRadLen(0,0);
204 if (!gGeoManager) {
205 gGeoManager = new TGeoManager("Geometry", "default geometry");
206 }
207 if (fZ - Int_t(fZ) > 1E-3)
208 Warning("ctor", "Material %s defined with fractional Z=%f", GetName(), fZ);
209 if (GetElement()) GetElement()->SetUsed();
211}
212
213////////////////////////////////////////////////////////////////////////////////
214/// Constructor.
215///
216/// \param name material name.
217/// \param elem
218/// \param rho material density in g/cm3.
219
221 :TNamed(name, ""), TAttFill(),
222 fIndex(0),
223 fA(0.),
224 fZ(0.),
225 fDensity(rho),
226 fRadLen(0.),
227 fIntLen(0.),
228 fTemperature(0.),
229 fPressure(0.),
230 fState(kMatStateUndefined),
231 fShader(NULL),
232 fCerenkov(NULL),
233 fElement(elem),
234 fUserExtension(0),
235 fFWExtension(0)
236{
237 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
238 fName = fName.Strip();
240 fIndex = -1;
241 fA = elem->A();
242 fZ = elem->Z();
243 SetRadLen(0,0);
247 if (!gGeoManager) {
248 gGeoManager = new TGeoManager("Geometry", "default geometry");
249 }
250 if (fZ - Int_t(fZ) > 1E-3)
251 Warning("ctor", "Material %s defined with fractional Z=%f", GetName(), fZ);
252 if (GetElement()) GetElement()->SetUsed();
254}
255
256////////////////////////////////////////////////////////////////////////////////
257
259 TNamed(gm),
260 TAttFill(gm),
261 fIndex(gm.fIndex),
262 fA(gm.fA),
263 fZ(gm.fZ),
264 fDensity(gm.fDensity),
265 fRadLen(gm.fRadLen),
266 fIntLen(gm.fIntLen),
267 fTemperature(gm.fTemperature),
268 fPressure(gm.fPressure),
269 fState(gm.fState),
270 fShader(gm.fShader),
271 fCerenkov(gm.fCerenkov),
272 fElement(gm.fElement),
273 fUserExtension(gm.fUserExtension->Grab()),
274 fFWExtension(gm.fFWExtension->Grab())
275
276{
277 //copy constructor
278 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
280 TIter next(&fProperties);
282 while ((property = (TNamed*)next())) fProperties.Add(new TNamed(*property));
283}
284
285////////////////////////////////////////////////////////////////////////////////
286///assignment operator
287
289{
290 if(this!=&gm) {
293 fIndex=gm.fIndex;
294 fA=gm.fA;
295 fZ=gm.fZ;
297 fRadLen=gm.fRadLen;
298 fIntLen=gm.fIntLen;
301 fState=gm.fState;
302 fShader=gm.fShader;
308 TIter next(&fProperties);
310 while ((property = (TNamed*)next())) fProperties.Add(new TNamed(*property));
311 }
312 return *this;
313}
314
315////////////////////////////////////////////////////////////////////////////////
316/// Destructor
317
319{
322}
323
324////////////////////////////////////////////////////////////////////////////////
325/// Connect user-defined extension to the material. The material "grabs" a copy, so
326/// the original object can be released by the producer. Release the previously
327/// connected extension if any.
328///
329/// NOTE: This interface is intended for user extensions and is guaranteed not
330/// to be used by TGeo
331
333{
335 fUserExtension = 0;
336 if (ext) fUserExtension = ext->Grab();
337}
338
339//_____________________________________________________________________________
340const char *TGeoMaterial::GetPropertyRef(const char *property) const
341{
342 // Find reference for a given property
344 return (prop) ? prop->GetTitle() : nullptr;
345}
346
347//_____________________________________________________________________________
349{
350 // Find reference for a given property
352 if ( !prop ) return nullptr;
353 return gGeoManager->GetGDMLMatrix(prop->GetTitle());
354}
355
356//_____________________________________________________________________________
358{
359 // Find reference for a given property
361 if ( !prop ) return nullptr;
362 return gGeoManager->GetGDMLMatrix(prop->GetTitle());
363}
364
365//_____________________________________________________________________________
366const char *TGeoMaterial::GetConstPropertyRef(const char *property) const
367{
368 // Find reference for a given constant property
370 return (prop) ? prop->GetTitle() : nullptr;
371}
372
373//_____________________________________________________________________________
375{
376 // Find reference for a given constant property
378 if (!prop) {
379 if (err) *err = kTRUE;
380 return 0.;
381 }
382 return gGeoManager->GetProperty(prop->GetTitle(), err);
383}
384
385//_____________________________________________________________________________
387{
388 // Find reference for a given constant property
390 if (!prop) {
391 if (err) *err = kTRUE;
392 return 0.;
393 }
394 return gGeoManager->GetProperty(prop->GetTitle(), err);
395}
396
397//_____________________________________________________________________________
398bool TGeoMaterial::AddProperty(const char *property, const char *ref)
399{
402 Error("AddProperty", "Property %s already added to material %s",
403 property, GetName());
404 return false;
405 }
406 fProperties.Add(new TNamed(property, ref));
407 return true;
408}
409
410//_____________________________________________________________________________
411bool TGeoMaterial::AddConstProperty(const char *property, const char *ref)
412{
415 Error("AddConstProperty", "Constant property %s already added to material %s",
416 property, GetName());
417 return false;
418 }
420 return true;
421}
422
423////////////////////////////////////////////////////////////////////////////////
424/// Connect framework defined extension to the material. The material "grabs" a copy,
425/// so the original object can be released by the producer. Release the previously
426/// connected extension if any.
427///
428/// NOTE: This interface is intended for the use by TGeo and the users should
429/// NOT connect extensions using this method
430
432{
434 fFWExtension = 0;
435 if (ext) fFWExtension = ext->Grab();
436}
437
438////////////////////////////////////////////////////////////////////////////////
439/// Get a copy of the user extension pointer. The user must call Release() on
440/// the copy pointer once this pointer is not needed anymore (equivalent to
441/// delete() after calling new())
442
444{
445 if (fUserExtension) return fUserExtension->Grab();
446 return 0;
447}
448
449////////////////////////////////////////////////////////////////////////////////
450/// Get a copy of the framework extension pointer. The user must call Release() on
451/// the copy pointer once this pointer is not needed anymore (equivalent to
452/// delete() after calling new())
453
455{
456 if (fFWExtension) return fFWExtension->Grab();
457 return 0;
458}
459
460////////////////////////////////////////////////////////////////////////////////
461/// Provide a pointer name containing uid.
462
464{
465 static TString name;
466 name = TString::Format("pMat%d", GetUniqueID());
467 return (char*)name.Data();
468}
469
470////////////////////////////////////////////////////////////////////////////////
471/// Set radiation/absorption lengths. If the values are negative, their absolute value
472/// is taken, otherwise radlen is recomputed using G3 formula.
473
475{
476 fRadLen = TMath::Abs(radlen);
477 fIntLen = TMath::Abs(intlen);
478 // Check for vacuum
479 if (fA<0.9 || fZ<0.9) {
480 if (radlen<-1e5 || intlen<-1e-5) {
481 Error("SetRadLen","Material %s: user values taken for vacuum: radlen=%g or intlen=%g - too small", GetName(),fRadLen, fIntLen);
482 return;
483 }
484 // Ignore positive values and take big numbers
485 if (radlen>=0) fRadLen = 1.E30;
486 if (intlen>=0) fIntLen = 1.E30;
487 return;
488 }
490 // compute radlen systematically with G3 formula for a valid material
491 if ( typ == TGeoManager::kRootUnits && radlen>=0 ) {
492 //taken grom Geant3 routine GSMATE
493 constexpr Double_t alr2av = 1.39621E-03*TGeoUnit::cm2;
494 constexpr Double_t al183 = 5.20948;
498 }
499 else if ( typ == TGeoManager::kG4Units && radlen>=0 ) {
500 //taken grom Geant3 routine GSMATE
501 constexpr Double_t alr2av = 1.39621E-03*TGeant4Unit::cm2;
502 constexpr Double_t al183 = 5.20948;
506 }
507 // Compute interaction length using the same formula as in GEANT4
508 if ( typ == TGeoManager::kRootUnits && intlen>=0 ) {
509 constexpr Double_t lambda0 = 35.*TGeoUnit::g/TGeoUnit::cm2; // [g/cm^2]
510 Double_t nilinv = 0.0;
511 TGeoElement *elem = GetElement();
512 if (!elem) {
513 Fatal("SetRadLen", "Element not found for material %s", GetName());
514 return;
515 }
516 Double_t nbAtomsPerVolume = TGeoUnit::Avogadro*fDensity/elem->A();
517 nilinv += nbAtomsPerVolume*TMath::Power(elem->Neff(), 0.6666667);
518 nilinv *= TGeoUnit::amu/lambda0;
519 fIntLen = (nilinv<=0) ? TGeoShape::Big() : (TGeoUnit::cm/nilinv);
520 }
521 else if ( typ == TGeoManager::kG4Units && intlen>=0 ) {
522 constexpr Double_t lambda0 = 35.*TGeant4Unit::g/TGeant4Unit::cm2; // [g/cm^2]
523 Double_t nilinv = 0.0;
524 TGeoElement *elem = GetElement();
525 if (!elem) {
526 Fatal("SetRadLen", "Element not found for material %s", GetName());
527 return;
528 }
529 Double_t nbAtomsPerVolume = TGeant4Unit::Avogadro*fDensity/elem->A();
530 nilinv += nbAtomsPerVolume*TMath::Power(elem->Neff(), 0.6666667);
531 nilinv *= TGeant4Unit::amu/lambda0;
532 fIntLen = (nilinv<=0) ? TGeoShape::Big() : (TGeant4Unit::cm/nilinv);
533 }
534}
535
536////////////////////////////////////////////////////////////////////////////////
537/// static function
538/// Compute Coulomb correction for pair production and Brem
539/// REFERENCE : EGS MANUAL SLAC 210 - UC32 - JUNE 78
540/// FORMULA 2.7.17
541
543{
546 Double_t az2 = az*az;
547 Double_t az4 = az2 * az2;
548 Double_t fp = ( 0.0083*az4 + 0.20206 + 1./(1.+az2) ) * az2;
549 Double_t fm = ( 0.0020*az4 + 0.0369 ) * az4;
550 return fp - fm;
551}
552
553////////////////////////////////////////////////////////////////////////////////
554/// return true if the other material has the same physical properties
555
557{
558 if (other==this) return kTRUE;
559 if (other->IsMixture()) return kFALSE;
560 if (TMath::Abs(fA-other->GetA())>1E-3) return kFALSE;
561 if (TMath::Abs(fZ-other->GetZ())>1E-3) return kFALSE;
562 if (TMath::Abs(fDensity-other->GetDensity())>1E-6) return kFALSE;
563 if (GetCerenkovProperties() != other->GetCerenkovProperties()) return kFALSE;
564// if (fRadLen != other->GetRadLen()) return kFALSE;
565// if (fIntLen != other->GetIntLen()) return kFALSE;
566 return kTRUE;
567}
568
569////////////////////////////////////////////////////////////////////////////////
570/// print characteristics of this material
571
572void TGeoMaterial::Print(const Option_t * /*option*/) const
573{
574 printf("Material %s %s A=%g Z=%g rho=%g radlen=%g intlen=%g index=%i\n", GetName(), GetTitle(),
576}
577
578////////////////////////////////////////////////////////////////////////////////
579/// Save a primitive as a C++ statement(s) on output stream "out".
580
581void TGeoMaterial::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
582{
584 char *name = GetPointerName();
585 out << "// Material: " << GetName() << std::endl;
586 out << " a = " << fA << ";" << std::endl;
587 out << " z = " << fZ << ";" << std::endl;
588 out << " density = " << fDensity << ";" << std::endl;
589 out << " radl = " << fRadLen << ";" << std::endl;
590 out << " absl = " << fIntLen << ";" << std::endl;
591
592 out << " " << name << " = new TGeoMaterial(\"" << GetName() << "\", a,z,density,radl,absl);" << std::endl;
593 out << " " << name << "->SetIndex(" << GetIndex() << ");" << std::endl;
595}
596
597////////////////////////////////////////////////////////////////////////////////
598/// Get some default color related to this material.
599
601{
603 return (2+id%6);
604}
605
606////////////////////////////////////////////////////////////////////////////////
607/// Get a pointer to the element this material is made of.
608/// This second call is to avoid warnings to not call a virtual
609/// method from the constructor
610
612{
613 if (fElement) return fElement;
615 return table->GetElement(Int_t(fZ));
616}
617
618////////////////////////////////////////////////////////////////////////////////
619/// Get a pointer to the element this material is made of.
620
622{
623 if (fElement) return fElement;
625 return table->GetElement(Int_t(fZ));
626}
627
628////////////////////////////////////////////////////////////////////////////////
629/// Single interface to get element properties.
630
632{
633 a = fA;
634 z = fZ;
635 w = 1.;
636}
637
638////////////////////////////////////////////////////////////////////////////////
639/// Retrieve material index in the list of materials
640
642{
643 if (fIndex>=0) return fIndex;
645 fIndex = matlist->IndexOf(this);
646 return fIndex;
647}
648
649////////////////////////////////////////////////////////////////////////////////
650/// Create the material representing the decay product of this material at a
651/// given time. The precision represent the minimum cumulative branching ratio for
652/// which decay products are still taken into account.
653
655{
656 TObjArray *pop = new TObjArray();
657 if (!fElement || !fElement->IsRadioNuclide()) return this;
658 FillMaterialEvolution(pop, precision);
659 Int_t ncomp = pop->GetEntriesFast();
660 if (!ncomp) return this;
661 TGeoElementRN *el;
662 Double_t *weight = new Double_t[ncomp];
663 Double_t amed = 0.;
664 Int_t i;
665 for (i=0; i<ncomp; i++) {
666 el = (TGeoElementRN *)pop->At(i);
667 weight[i] = el->Ratio()->Concentration(time) * el->A();
668 amed += weight[i];
669 }
670 Double_t rho = fDensity*amed/fA;
671 TGeoMixture *mix = 0;
672 Int_t ncomp1 = ncomp;
673 for (i=0; i<ncomp; i++) {
674 if ((weight[i]/amed)<precision) {
675 amed -= weight[i];
676 ncomp1--;
677 }
678 }
679 if (ncomp1<2) {
680 el = (TGeoElementRN *)pop->At(0);
681 delete [] weight;
682 delete pop;
683 if (ncomp1==1) return new TGeoMaterial(TString::Format("%s-evol",GetName()), el, rho);
684 return NULL;
685 }
686 mix = new TGeoMixture(TString::Format("%s-evol",GetName()), ncomp, rho);
687 for (i=0; i<ncomp; i++) {
688 weight[i] /= amed;
689 if (weight[i]<precision) continue;
690 el = (TGeoElementRN *)pop->At(i);
691 mix->AddElement(el, weight[i]);
692 }
693 delete [] weight;
694 delete pop;
695 return mix;
696}
697
698////////////////////////////////////////////////////////////////////////////////
699/// Fills a user array with all the elements deriving from the possible
700/// decay of the top element composing the mixture. Each element contained
701/// by `<population>` may be a radionuclide having a Bateman solution attached.
702/// The precision represent the minimum cumulative branching ratio for
703/// which decay products are still taken into account.
704/// To visualize the time evolution of each decay product one can use:
705/// ~~~ {.cpp}
706/// TGeoElement *elem = population->At(index);
707/// TGeoElementRN *elemrn = 0;
708/// if (elem->IsRadioNuclide()) elemrn = (TGeoElementRN*)elem;
709/// ~~~
710/// One can get Ni/N1(t=0) at any moment of time. Ni is the number of atoms
711/// of one of the decay products, N1(0) is the number of atoms of the top
712/// element at t=0.
713/// ~~~ {.cpp}
714/// Double_t fraction_weight = elemrn->Ratio()->Concentration(time);
715/// ~~~
716/// One can also display the time evolution of the fractional weight:
717/// ~~~ {.cpp}
718/// elemrn->Ratio()->Draw(option);
719/// ~~~
720
722{
723 if (population->GetEntriesFast()) {
724 Error("FillMaterialEvolution", "Provide an empty array !");
725 return;
726 }
728 TGeoElement *elem;
729 TGeoElementRN *elemrn;
730 TIter next(table->GetElementsRN());
731 while ((elemrn=(TGeoElementRN*)next())) elemrn->ResetRatio();
732 elem = GetElement();
733 if (!elem) {
734 Fatal("FillMaterialEvolution", "Element not found for material %s", GetName());
735 return;
736 }
737 if (!elem->IsRadioNuclide()) {
738 population->Add(elem);
739 return;
740 }
741 elemrn = (TGeoElementRN*)elem;
742 elemrn->FillPopulation(population, precision);
743}
744
745/** \class TGeoMixture
746\ingroup Materials_classes
747
748Mixtures of elements.
749
750*/
751
753
754////////////////////////////////////////////////////////////////////////////////
755/// Default constructor
756
758{
759 fNelements = 0;
760 fZmixture = 0;
761 fAmixture = 0;
762 fWeights = 0;
763 fNatoms = 0;
765 fElements = 0;
766}
767
768////////////////////////////////////////////////////////////////////////////////
769/// constructor
770
771TGeoMixture::TGeoMixture(const char *name, Int_t /*nel*/, Double_t rho)
773{
774 fZmixture = 0;
775 fAmixture = 0;
776 fWeights = 0;
777 fNelements = 0;
778 fNatoms = 0;
780 fDensity = rho;
781 fElements = 0;
782 if (fDensity < 0) fDensity = 0.001;
783}
784
785////////////////////////////////////////////////////////////////////////////////
786/// Destructor
787
789{
790 if (fZmixture) delete[] fZmixture;
791 if (fAmixture) delete[] fAmixture;
792 if (fWeights) delete[] fWeights;
793 if (fNatoms) delete[] fNatoms;
795 if (fElements) delete fElements;
796}
797
798////////////////////////////////////////////////////////////////////////////////
799/// Compute effective A/Z and radiation length
800
802{
806 const Double_t amu = (typ==TGeoManager::kRootUnits) ? TGeoUnit::amu : TGeant4Unit::amu; // [MeV/c^2]
809 const Double_t alr2av = 1.39621E-03 * cm2;
810 const Double_t al183 = 5.20948;
811 const Double_t lambda0 = 35.*gram/cm2; // [g/cm^2]
812 Double_t radinv = 0.0;
813 Double_t nilinv = 0.0;
814 Double_t nbAtomsPerVolume;
815 fA = 0;
816 fZ = 0;
817 for (Int_t j=0;j<fNelements;j++) {
818 if (fWeights[j] <= 0) continue;
819 fA += fWeights[j]*fAmixture[j];
820 fZ += fWeights[j]*fZmixture[j];
821 nbAtomsPerVolume = na*fDensity*fWeights[j]/GetElement(j)->A();
822 nilinv += nbAtomsPerVolume*TMath::Power(GetElement(j)->Neff(), 0.6666667);
823 Double_t zc = fZmixture[j];
824 Double_t alz = TMath::Log(zc)/3.;
825 Double_t xinv = zc*(zc+TGeoMaterial::ScreenFactor(zc))*
826 (al183-alz-TGeoMaterial::Coulomb(zc))/fAmixture[j];
827 radinv += xinv*fWeights[j];
828 }
829 radinv *= alr2av*fDensity;
830 if (radinv > 0) fRadLen = cm/radinv;
831 // Compute interaction length
832 nilinv *= amu/lambda0;
833 fIntLen = (nilinv<=0) ? TGeoShape::Big() : (cm/nilinv);
834}
835
836////////////////////////////////////////////////////////////////////////////////
837/// add an element to the mixture using fraction by weight
838/// Check if the element is already defined
839
841{
843
844 // Check preconditions
845 if (weight < 0e0) {
846 Fatal("AddElement", "Cannot add element with negative weight %g to mixture %s", weight, GetName());
847 }
848 else if ( weight < std::numeric_limits<Double_t>::epsilon() ) {
849 return;
850 }
851 else if (z<1 || z>table->GetNelements()-1) {
852 Fatal("AddElement", "Cannot add element having Z=%d to mixture %s", (Int_t)z, GetName());
853 }
854 Int_t i;
855 for (i=0; i<fNelements; i++) {
856 if (!fElements && TMath::Abs(z-fZmixture[i])<1.e-6 && TMath::Abs(a-fAmixture[i])<1.e-6) {
857 fWeights[i] += weight;
859 return;
860 }
861 }
862 if (!fNelements) {
863 fZmixture = new Double_t[1];
864 fAmixture = new Double_t[1];
865 fWeights = new Double_t[1];
866 } else {
867 Int_t nelements = fNelements+1;
868 Double_t *zmixture = new Double_t[nelements];
869 Double_t *amixture = new Double_t[nelements];
870 Double_t *weights = new Double_t[nelements];
871 for (Int_t j=0; j<fNelements; j++) {
872 zmixture[j] = fZmixture[j];
873 amixture[j] = fAmixture[j];
874 weights[j] = fWeights[j];
875 }
876 delete [] fZmixture;
877 delete [] fAmixture;
878 delete [] fWeights;
879 fZmixture = zmixture;
880 fAmixture = amixture;
881 fWeights = weights;
882 }
883
884 fNelements++;
885 i = fNelements - 1;
886 fZmixture[i] = z;
887 fAmixture[i] = a;
888 fWeights[i] = weight;
889 if (z - Int_t(z) > 1E-3)
890 Warning("DefineElement", "Mixture %s has element defined with fractional Z=%f", GetName(), z);
892 table->GetElement((Int_t)z)->SetDefined();
893
894 //compute equivalent radiation length (taken from Geant3/GSMIXT)
896}
897
898////////////////////////////////////////////////////////////////////////////////
899/// Define one component of the mixture as an existing material/mixture.
900
902{
903 TGeoElement *elnew, *elem;
904 Double_t a,z;
905
906 // Check preconditions
907 if (!mat) {
908 Fatal("AddElement", "Cannot add INVALID material to mixture %s", GetName());
909 }
910 else if (weight < 0e0) {
911 Fatal("AddElement", "Cannot add material %s with negative weight %g to mixture %s",
912 mat->GetName(), weight, GetName());
913 }
914 else if ( weight < std::numeric_limits<Double_t>::epsilon() ) {
915 return;
916 }
917 if (!mat->IsMixture()) {
918 elem = mat->GetBaseElement();
919 if (elem) {
920 AddElement(elem, weight);
921 } else {
922 a = mat->GetA();
923 z = mat->GetZ();
924 AddElement(a, z, weight);
925 }
926 return;
927 }
928 // The material is a mixture.
929 TGeoMixture *mix = (TGeoMixture*)mat;
930 Double_t wnew;
931 Int_t nelem = mix->GetNelements();
932 Bool_t elfound;
933 Int_t i,j;
934 // loop the elements of the daughter mixture
935 for (i=0; i<nelem; i++) {
936 elfound = kFALSE;
937 elnew = mix->GetElement(i);
938 if (!elnew) continue;
939 // check if we have the element already defined in the parent mixture
940 for (j=0; j<fNelements; j++) {
941 if ( fWeights[j] < 0e0 ) continue;
942 elem = GetElement(j);
943 if (elem == elnew) {
944 // element found, compute new weight
945 fWeights[j] += weight * (mix->GetWmixt())[i];
946 elfound = kTRUE;
947 break;
948 }
949 }
950 if (elfound) continue;
951 // element not found, define it
952 wnew = weight * (mix->GetWmixt())[i];
953 AddElement(elnew, wnew);
954 }
955}
956
957////////////////////////////////////////////////////////////////////////////////
958/// add an element to the mixture using fraction by weight
959
961{
962 TGeoElement *elemold;
964 if (!fElements) fElements = new TObjArray(128);
965 Bool_t exist = kFALSE;
966
967 // Check preconditions
968 if (!elem) {
969 Fatal("AddElement", "Cannot add INVALID element to mixture %s", GetName());
970 }
971 else if (weight < 0e0) {
972 Fatal("AddElement", "Cannot add element %s with negative weight %g to mixture %s",
973 elem->GetName(), weight, GetName());
974 }
975 else if ( weight < std::numeric_limits<Double_t>::epsilon() ) {
976 return;
977 }
978 // If previous elements were defined by A/Z, add corresponding TGeoElements
979 for (Int_t i=0; i<fNelements; i++) {
980 elemold = (TGeoElement*)fElements->At(i);
981 if (!elemold) {
982 fElements->AddAt(elemold = table->GetElement((Int_t)fZmixture[i]), i);
983 }
984 if (elemold == elem) exist = kTRUE;
985 }
986 if (!exist) {
988 }
989 AddElement(elem->A(), elem->Z(), weight);
990}
991
992////////////////////////////////////////////////////////////////////////////////
993/// Add a mixture element by number of atoms in the chemical formula.
994
996{
997 Int_t i,j;
998 Double_t amol;
999 TGeoElement *elemold;
1001 if (!fElements) fElements = new TObjArray(128);
1002 // Check if the element is already defined
1003 for (i=0; i<fNelements; i++) {
1004 elemold = (TGeoElement*)fElements->At(i);
1005 if (!elemold) fElements->AddAt(table->GetElement((Int_t)fZmixture[i]), i);
1006 else if (elemold != elem) continue;
1007 if ((elem==elemold) ||
1008 (TMath::Abs(elem->Z()-fZmixture[i])<1.e-6 && TMath::Abs(elem->A()-fAmixture[i])<1.e-6)) {
1009 fNatoms[i] += natoms;
1010 amol = 0.;
1011 for (j=0; j<fNelements; j++) amol += fAmixture[j]*fNatoms[j];
1012 for (j=0; j<fNelements; j++) fWeights[j] = fNatoms[j]*fAmixture[j]/amol;
1014 return;
1015 }
1016 }
1017 // New element
1018 if (!fNelements) {
1019 fZmixture = new Double_t[1];
1020 fAmixture = new Double_t[1];
1021 fWeights = new Double_t[1];
1022 fNatoms = new Int_t[1];
1023 } else {
1024 if (!fNatoms) {
1025 Fatal("AddElement", "Cannot add element by natoms in mixture %s after defining elements by weight",
1026 GetName());
1027 return;
1028 }
1029 Int_t nelements = fNelements+1;
1030 Double_t *zmixture = new Double_t[nelements];
1031 Double_t *amixture = new Double_t[nelements];
1032 Double_t *weights = new Double_t[nelements];
1033 Int_t *nnatoms = new Int_t[nelements];
1034 for (j=0; j<fNelements; j++) {
1035 zmixture[j] = fZmixture[j];
1036 amixture[j] = fAmixture[j];
1037 weights[j] = fWeights[j];
1038 nnatoms[j] = fNatoms[j];
1039 }
1040 delete [] fZmixture;
1041 delete [] fAmixture;
1042 delete [] fWeights;
1043 delete [] fNatoms;
1044 fZmixture = zmixture;
1045 fAmixture = amixture;
1046 fWeights = weights;
1047 fNatoms = nnatoms;
1048 }
1049 fNelements++;
1050 Int_t iel = fNelements-1;
1051 fZmixture[iel] = elem->Z();
1052 fAmixture[iel] = elem->A();
1053 fNatoms[iel] = natoms;
1054 fElements->AddAtAndExpand(elem, iel);
1055 amol = 0.;
1056 for (i=0; i<fNelements; i++) {
1057 if (fNatoms[i]<=0) return;
1058 amol += fAmixture[i]*fNatoms[i];
1059 }
1060 for (i=0; i<fNelements; i++) fWeights[i] = fNatoms[i]*fAmixture[i]/amol;
1061 table->GetElement(elem->Z())->SetDefined();
1063}
1064
1065////////////////////////////////////////////////////////////////////////////////
1066/// Define the mixture element at index iel by number of atoms in the chemical formula.
1067
1069{
1071 TGeoElement *elem = table->GetElement(z);
1072 if (!elem) {
1073 Fatal("DefineElement", "In mixture %s, element with Z=%i not found",GetName(),z);
1074 return;
1075 }
1076 AddElement(elem, natoms);
1077}
1078
1079////////////////////////////////////////////////////////////////////////////////
1080/// Retrieve the pointer to the element corresponding to component I.
1081
1083{
1084 if (i<0 || i>=fNelements) {
1085 Error("GetElement", "Mixture %s has only %d elements", GetName(), fNelements);
1086 return 0;
1087 }
1088 TGeoElement *elem = 0;
1089 if (fElements) elem = (TGeoElement*)fElements->At(i);
1090 if (elem) return elem;
1092 return table->GetElement(Int_t(fZmixture[i]));
1093}
1094
1095////////////////////////////////////////////////////////////////////////////////
1096/// Get specific activity (in Bq/gram) for the whole mixture (no argument) or
1097/// for a given component.
1098
1100{
1101 if (i>=0 && i<fNelements) return fWeights[i]*GetElement(i)->GetSpecificActivity();
1102 Double_t sa = 0;
1103 for (Int_t iel=0; iel<fNelements; iel++) {
1104 sa += fWeights[iel]*GetElement(iel)->GetSpecificActivity();
1105 }
1106 return sa;
1107}
1108
1109////////////////////////////////////////////////////////////////////////////////
1110/// Return true if the other material has the same physical properties
1111
1113{
1114 if (other->IsEqual(this)) return kTRUE;
1115 if (!other->IsMixture()) return kFALSE;
1116 TGeoMixture *mix = (TGeoMixture*)other;
1117 if (!mix) return kFALSE;
1118 if (fNelements != mix->GetNelements()) return kFALSE;
1119 if (TMath::Abs(fA-other->GetA())>1E-3) return kFALSE;
1120 if (TMath::Abs(fZ-other->GetZ())>1E-3) return kFALSE;
1121 if (TMath::Abs(fDensity-other->GetDensity())>1E-6) return kFALSE;
1122 if (GetCerenkovProperties() != other->GetCerenkovProperties()) return kFALSE;
1123// if (fRadLen != other->GetRadLen()) return kFALSE;
1124// if (fIntLen != other->GetIntLen()) return kFALSE;
1125 for (Int_t i=0; i<fNelements; i++) {
1126 if (TMath::Abs(fZmixture[i]-(mix->GetZmixt())[i])>1E-3) return kFALSE;
1127 if (TMath::Abs(fAmixture[i]-(mix->GetAmixt())[i])>1E-3) return kFALSE;
1128 if (TMath::Abs(fWeights[i]-(mix->GetWmixt())[i])>1E-3) return kFALSE;
1129 }
1130 return kTRUE;
1131}
1132
1133////////////////////////////////////////////////////////////////////////////////
1134/// print characteristics of this material
1135
1136void TGeoMixture::Print(const Option_t * /*option*/) const
1137{
1138 printf("Mixture %s %s Aeff=%g Zeff=%g rho=%g radlen=%g intlen=%g index=%i\n", GetName(), GetTitle(),
1140 for (Int_t i=0; i<fNelements; i++) {
1141 if (fElements && fElements->At(i)) {
1142 fElements->At(i)->Print();
1143 continue;
1144 }
1145 if (fNatoms) printf(" Element #%i : %s Z=%6.2f A=%6.2f w=%6.3f natoms=%d\n", i, GetElement(i)->GetName(),fZmixture[i],
1146 fAmixture[i], fWeights[i], fNatoms[i]);
1147 else printf(" Element #%i : %s Z=%6.2f A=%6.2f w=%6.3f\n", i, GetElement(i)->GetName(),fZmixture[i],
1148 fAmixture[i], fWeights[i]);
1149 }
1150}
1151
1152////////////////////////////////////////////////////////////////////////////////
1153/// Save a primitive as a C++ statement(s) on output stream "out".
1154
1155void TGeoMixture::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
1156{
1158 char *name = GetPointerName();
1159 out << "// Mixture: " << GetName() << std::endl;
1160 out << " nel = " << fNelements << ";" << std::endl;
1161 out << " density = " << fDensity << ";" << std::endl;
1162 out << " " << name << " = new TGeoMixture(\"" << GetName() << "\", nel,density);" << std::endl;
1163 for (Int_t i=0; i<fNelements; i++) {
1164 TGeoElement *el = GetElement(i);
1165 out << " a = " << fAmixture[i] << "; z = "<< fZmixture[i] << "; w = " << fWeights[i] << "; // " << el->GetName() << std::endl;
1166 out << " " << name << "->DefineElement(" << i << ",a,z,w);" << std::endl;
1167 }
1168 out << " " << name << "->SetIndex(" << GetIndex() << ");" << std::endl;
1170}
1171
1172////////////////////////////////////////////////////////////////////////////////
1173/// Create the mixture representing the decay product of this material at a
1174/// given time. The precision represent the minimum cumulative branching ratio for
1175/// which decay products are still taken into account.
1176
1178{
1179 TObjArray *pop = new TObjArray();
1180 FillMaterialEvolution(pop, precision);
1181 Int_t ncomp = pop->GetEntriesFast();
1182 if (!ncomp) return this;
1183 TGeoElement *elem;
1184 TGeoElementRN *el;
1185 Double_t *weight = new Double_t[ncomp];
1186 Double_t amed = 0.;
1187 Int_t i, j;
1188 for (i=0; i<ncomp; i++) {
1189 elem = (TGeoElement *)pop->At(i);
1190 if (!elem->IsRadioNuclide()) {
1191 j = fElements->IndexOf(elem);
1192 weight[i] = fWeights[j]*fAmixture[0]/fWeights[0];
1193 } else {
1194 el = (TGeoElementRN*)elem;
1195 weight[i] = el->Ratio()->Concentration(time) * el->A();
1196 }
1197 amed += weight[i];
1198 }
1199 Double_t rho = fDensity * fWeights[0] * amed/fAmixture[0];
1200 TGeoMixture *mix = 0;
1201 Int_t ncomp1 = ncomp;
1202 for (i=0; i<ncomp; i++) {
1203 if ((weight[i]/amed)<precision) {
1204 amed -= weight[i];
1205 ncomp1--;
1206 }
1207 }
1208 if (ncomp1<2) {
1209 el = (TGeoElementRN *)pop->At(0);
1210 delete [] weight;
1211 delete pop;
1212 if (ncomp1==1) return new TGeoMaterial(TString::Format("%s-evol",GetName()), el, rho);
1213 return NULL;
1214 }
1215 mix = new TGeoMixture(TString::Format("%s-evol",GetName()), ncomp, rho);
1216 for (i=0; i<ncomp; i++) {
1217 weight[i] /= amed;
1218 if (weight[i]<precision) continue;
1219 el = (TGeoElementRN *)pop->At(i);
1220 mix->AddElement(el, weight[i]);
1221 }
1222 delete [] weight;
1223 delete pop;
1224 return mix;
1225}
1226
1227////////////////////////////////////////////////////////////////////////////////
1228/// Fills a user array with all the elements deriving from the possible
1229/// decay of the top elements composing the mixture. Each element contained
1230/// by `<population>` may be a radionuclide having a Bateman solution attached.
1231/// The precision represent the minimum cumulative branching ratio for
1232/// which decay products are still taken into account.
1233/// To visualize the time evolution of each decay product one can use:
1234/// ~~~ {.cpp}
1235/// TGeoElement *elem = population->At(index);
1236/// TGeoElementRN *elemrn = 0;
1237/// if (elem->IsRadioNuclide()) elemrn = (TGeoElementRN*)elem;
1238/// ~~~
1239/// One can get Ni/N1(t=0) at any moment of time. Ni is the number of atoms
1240/// of one of the decay products, N1(0) is the number of atoms of the first top
1241/// element at t=0.
1242/// ~~~ {.cpp}
1243/// Double_t fraction_weight = elemrn->Ratio()->Concentration(time);
1244/// ~~~
1245/// One can also display the time evolution of the fractional weight:
1246/// ~~~ {.cpp}
1247/// elemrn->Ratio()->Draw(option);
1248/// ~~~
1249
1251{
1252 if (population->GetEntriesFast()) {
1253 Error("FillMaterialEvolution", "Provide an empty array !");
1254 return;
1255 }
1257 TGeoElement *elem;
1258 TGeoElementRN *elemrn;
1259 TIter next(table->GetElementsRN());
1260 while ((elemrn=(TGeoElementRN*)next())) elemrn->ResetRatio();
1261 Double_t factor;
1262 for (Int_t i=0; i<fNelements; i++) {
1263 elem = GetElement(i);
1264 if (!elem->IsRadioNuclide()) {
1265 population->Add(elem);
1266 continue;
1267 }
1268 elemrn = (TGeoElementRN*)elem;
1269 factor = fWeights[i]*fAmixture[0]/(fWeights[0]*fAmixture[i]);
1270 elemrn->FillPopulation(population, precision, factor);
1271 }
1272}
1273
1274////////////////////////////////////////////////////////////////////////////////
1275/// static function
1276/// Compute screening factor for pair production and Bremsstrahlung
1277/// REFERENCE : EGS MANUAL SLAC 210 - UC32 - JUNE 78
1278/// FORMULA 2.7.22
1279
1281{
1282 const Double_t al183= 5.20948 , al1440 = 7.27239;
1283 Double_t alz = TMath::Log(z)/3.;
1284 Double_t factor = (al1440 - 2*alz) / (al183 - alz - TGeoMaterial::Coulomb(z));
1285 return factor;
1286}
1287
1288////////////////////////////////////////////////////////////////////////////////
1289/// Compute Derived Quantities as in Geant4
1290
1292{
1295
1297
1299
1300 // Formula taken from G4Material.cxx L312
1301 for (Int_t i=0; i<fNelements; ++i) {
1303 }
1306}
1307
1308
1309////////////////////////////////////////////////////////////////////////////////
1310/// Compute Radiation Length based on Geant4 formula
1311
1313{
1314 // Formula taken from G4Material.cxx L556
1316 Double_t radinv = 0.0 ;
1317 for (Int_t i=0;i<fNelements;++i) {
1318 radinv += fVecNbOfAtomsPerVolume[i]*((TGeoElement*)fElements->At(i))->GetfRadTsai();
1319 }
1320 fRadLen = (radinv <= 0.0 ? DBL_MAX : cm/radinv);
1321}
1322
1323////////////////////////////////////////////////////////////////////////////////
1324/// Compute Nuclear Interaction Length based on Geant4 formula
1326{
1327 // Formula taken from G4Material.cxx L567
1332 const Double_t lambda0 = 35*g/(cm*cm);
1333 const Double_t twothird = 2.0/3.0;
1334 Double_t NILinv = 0.0;
1335 for (Int_t i=0; i<fNelements; ++i) {
1336 Int_t Z = static_cast<Int_t>(((TGeoElement*)fElements->At(i))->Z()+0.5);
1337 Double_t A = ((TGeoElement*)fElements->At(i))->Neff();
1338 if(1 == Z) {
1339 NILinv += fVecNbOfAtomsPerVolume[i]*A;
1340 } else {
1341 NILinv += fVecNbOfAtomsPerVolume[i]*TMath::Exp(twothird*TMath::Log(A));
1342 }
1343 }
1344 NILinv *= amu/lambda0;
1345 fIntLen = (NILinv <= 0.0 ? DBL_MAX : cm/NILinv);
1346}
#define e(i)
Definition: RSha256.hxx:103
int Int_t
Definition: RtypesCore.h:45
const Bool_t kFALSE
Definition: RtypesCore.h:101
const Bool_t kTRUE
Definition: RtypesCore.h:100
const char Option_t
Definition: RtypesCore.h:66
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 g
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
Definition: TGeoManager.h:602
ClassImp(TGeoMaterial)
static const Double_t STP_temperature
Definition: TGeoMaterial.h:31
static const Double_t STP_pressure
Definition: TGeoMaterial.h:32
Binding & operator=(OUT(*fun)(void))
Fill Area Attributes class.
Definition: TAttFill.h:19
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:34
Double_t Concentration(Double_t time) const
Find concentration of the element at a given time.
Class representing a radionuclidevoid TGeoManager::SetDefaultRootUnits() { if ( fgDefaultUnits == kRo...
Definition: TGeoElement.h:139
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.
TGeoBatemanSol * Ratio() const
Definition: TGeoElement.h:197
void ResetRatio()
Clears the existing ratio.
Table of elements.
Definition: TGeoElement.h:370
TGeoElement * GetElement(Int_t z)
Definition: TGeoElement.h:410
TObjArray * GetElementsRN() const
Definition: TGeoElement.h:413
Int_t GetNelements() const
Definition: TGeoElement.h:417
Base class for chemical elements.
Definition: TGeoElement.h:37
virtual Double_t GetSpecificActivity() const
Definition: TGeoElement.h:84
Double_t A() const
Definition: TGeoElement.h:76
void SetDefined(Bool_t flag=kTRUE)
Definition: TGeoElement.h:90
virtual Bool_t IsRadioNuclide() const
Definition: TGeoElement.h:87
Double_t Neff() const
Returns effective number of nucleons.
Int_t Z() const
Definition: TGeoElement.h:73
void SetUsed(Bool_t flag=kTRUE)
Definition: TGeoElement.h:91
ABC for user objects attached to TGeoVolume or TGeoNode.
Definition: TGeoExtension.h:20
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
Definition: TGeoManager.h:491
Base class describing materials.
Definition: TGeoMaterial.h:36
virtual ~TGeoMaterial()
Destructor.
Double_t GetConstProperty(const char *property, Bool_t *error=nullptr) const
void SetUserExtension(TGeoExtension *ext)
Connect user-defined extension to the material.
char * GetPointerName() const
Provide a pointer name containing uid.
virtual TObject * GetCerenkovProperties() const
Definition: TGeoMaterial.h:118
EGeoMaterialState fState
Definition: TGeoMaterial.h:58
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
virtual Bool_t IsMixture() const
Definition: TGeoMaterial.h:130
Double_t fPressure
Definition: TGeoMaterial.h:57
void SetFWExtension(TGeoExtension *ext)
Connect framework defined extension to the material.
bool AddConstProperty(const char *property, const char *ref)
Double_t fTemperature
Definition: TGeoMaterial.h:56
virtual void GetElementProp(Double_t &a, Double_t &z, Double_t &w, Int_t i=0)
Single interface to get element properties.
virtual void Print(const Option_t *option="") const
print characteristics of this material
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
Definition: TGeoMaterial.h:62
bool AddProperty(const char *property, const char *ref)
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Double_t fZ
Definition: TGeoMaterial.h:52
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
TGeoElement * GetBaseElement() const
Definition: TGeoMaterial.h:113
const char * GetPropertyRef(const char *property) const
Double_t fA
Definition: TGeoMaterial.h:51
TGDMLMatrix * GetProperty(const char *name) const
TObject * fCerenkov
Definition: TGeoMaterial.h:60
Double_t fDensity
Definition: TGeoMaterial.h:53
void SetUsed(Bool_t flag=kTRUE)
Definition: TGeoMaterial.h:139
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
Definition: TGeoMaterial.h:55
TGeoExtension * fUserExtension
Definition: TGeoMaterial.h:64
TList fConstProperties
Definition: TGeoMaterial.h:63
TGeoElement * fElement
Definition: TGeoMaterial.h:61
TGeoMaterial & operator=(const TGeoMaterial &)
assignment operator
TObject * fShader
Definition: TGeoMaterial.h:59
TGeoExtension * GrabUserExtension() const
Get a copy of the user extension pointer.
Double_t fRadLen
Definition: TGeoMaterial.h:54
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...
virtual Double_t GetA() const
Definition: TGeoMaterial.h:105
TGeoExtension * fFWExtension
Transient user-defined extension to materials.
Definition: TGeoMaterial.h:65
virtual Double_t GetDensity() const
Definition: TGeoMaterial.h:108
virtual Double_t GetZ() const
Definition: TGeoMaterial.h:106
Mixtures of elements.
Definition: TGeoMaterial.h:157
TObjArray * fElements
Definition: TGeoMaterial.h:166
virtual TGeoElement * GetElement(Int_t i=0) const
Retrieve the pointer to the element corresponding to component I.
void ComputeNuclearInterLength()
Compute Nuclear Interaction Length based on Geant4 formula.
Double_t * GetZmixt() const
Definition: TGeoMaterial.h:195
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
Definition: TGeoMaterial.h:161
virtual Double_t GetSpecificActivity(Int_t i=-1) const
Get specific activity (in Bq/gram) for the whole mixture (no argument) or for a given component.
TGeoMixture()
Default constructor.
virtual ~TGeoMixture()
Destructor.
void ComputeDerivedQuantities()
Compute Derived Quantities as in Geant4.
virtual void Print(const Option_t *option="") const
print characteristics of this material
void ComputeRadiationLength()
Compute Radiation Length based on Geant4 formula.
Double_t * fVecNbOfAtomsPerVolume
Definition: TGeoMaterial.h:165
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 elements composi...
Double_t * fAmixture
Definition: TGeoMaterial.h:162
void AverageProperties()
Compute effective A/Z and radiation length.
Int_t fNelements
Definition: TGeoMaterial.h:160
Double_t * GetWmixt() const
Definition: TGeoMaterial.h:197
virtual TGeoMaterial * DecayMaterial(Double_t time, Double_t precision=0.001)
Create the mixture representing the decay product of this material at a given time.
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
virtual Int_t GetNelements() const
Definition: TGeoMaterial.h:194
Int_t * fNatoms
Definition: TGeoMaterial.h:164
virtual Bool_t IsEq(const TGeoMaterial *other) const
Return true if the other material has the same physical properties.
Double_t * fWeights
Definition: TGeoMaterial.h:163
void DefineElement(Int_t iel, Double_t a, Double_t z, Double_t weight)
Definition: TGeoMaterial.h:215
Double_t * GetAmixt() const
Definition: TGeoMaterial.h:196
static Double_t Big()
Definition: TGeoShape.h:88
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:578
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:357
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:47
const char * GetTitle() const override
Returns title of object.
Definition: TNamed.h:48
TNamed()
Definition: TNamed.h:36
TString fName
Definition: TNamed.h:32
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition: TNamed.cxx:51
An array of TObjects.
Definition: TObjArray.h:31
Int_t GetEntriesFast() const
Definition: TObjArray.h:58
Int_t IndexOf(const TObject *obj) const override
Definition: TObjArray.cxx:605
void AddAt(TObject *obj, Int_t idx) override
Add object at position ids.
Definition: TObjArray.cxx:254
virtual void AddAtAndExpand(TObject *obj, Int_t idx)
Add object at position idx.
Definition: TObjArray.cxx:235
TObject * At(Int_t idx) const override
Definition: TObjArray.h:164
void Add(TObject *obj) override
Definition: TObjArray.h:68
virtual Bool_t IsEqual(const TObject *obj) const
Default equal comparison (objects are equal if they have the same address in memory).
Definition: TObject.cxx:485
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
Definition: TObject.h:187
virtual UInt_t GetUniqueID() const
Return the unique object id.
Definition: TObject.cxx:377
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:879
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
Definition: TObject.cxx:696
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:893
virtual void Fatal(const char *method, const char *msgfmt,...) const
Issue fatal error message.
Definition: TObject.cxx:921
virtual void Print(Option_t *option="") const
This method must be overridden when a class wants to print itself.
Definition: TObject.cxx:552
virtual Int_t IndexOf(const TObject *obj) const
Return index of object in collection.
Basic string class.
Definition: TString.h:136
TSubString Strip(EStripType s=kTrailing, char c=' ') const
Return a substring of self stripped at beginning and/or end.
Definition: TString.cxx:1131
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:2336
static double A[]
static constexpr double amu
static constexpr double fine_structure_const
static constexpr double gram
static constexpr double cm2
static constexpr double Avogadro
static constexpr double g
static constexpr double cm
static constexpr double amu
static constexpr double cm
static constexpr double fine_structure_const
static constexpr double gram
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:706
constexpr Double_t E()
Base of natural log: .
Definition: TMath.h:93
Double_t Log(Double_t x)
Returns the natural logarithm of x.
Definition: TMath.h:753
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Returns x raised to the power y.
Definition: TMath.h:718
constexpr Double_t Na()
Avogadro constant (Avogadro's Number) in .
Definition: TMath.h:284
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition: TMathBase.h:123
auto * a
Definition: textangle.C:12
double epsilon
Definition: triangle.c:618