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 Geometry_classes
14
15Base class describing materials.
16
17\image html geom_material.jpg
18*/
19
20#include <iostream>
21#include "TMath.h"
22#include "TObjArray.h"
23#include "TGeoElement.h"
24#include "TGeoManager.h"
25#include "TGeoExtension.h"
26#include "TGeoMaterial.h"
29#include "TGDMLMatrix.h"
30
31// statics and globals
32
34
35////////////////////////////////////////////////////////////////////////////////
36/// Default constructor
37
39 :TNamed(), TAttFill(),
40 fIndex(0),
41 fA(0.),
42 fZ(0.),
43 fDensity(0.),
44 fRadLen(0.),
45 fIntLen(0.),
46 fTemperature(0.),
47 fPressure(0.),
48 fState(kMatStateUndefined),
49 fShader(NULL),
50 fCerenkov(NULL),
51 fElement(NULL),
52 fUserExtension(0),
53 fFWExtension(0)
54{
55 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
57 fIndex = -1;
61}
62
63////////////////////////////////////////////////////////////////////////////////
64/// constructor
65
67 :TNamed(name, ""), TAttFill(),
68 fIndex(0),
69 fA(0.),
70 fZ(0.),
71 fDensity(0.),
72 fRadLen(0.),
73 fIntLen(0.),
74 fTemperature(0.),
75 fPressure(0.),
76 fState(kMatStateUndefined),
77 fShader(NULL),
78 fCerenkov(NULL),
79 fElement(NULL),
80 fUserExtension(0),
81 fFWExtension(0)
82{
83 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
84 fName = fName.Strip();
86 fIndex = -1;
90
91 if (!gGeoManager) {
92 gGeoManager = new TGeoManager("Geometry", "default geometry");
93 }
95}
96
97////////////////////////////////////////////////////////////////////////////////
98/// constructor
99
101 Double_t rho, Double_t radlen, Double_t intlen)
102 :TNamed(name, ""), TAttFill(),
103 fIndex(0),
104 fA(a),
105 fZ(z),
106 fDensity(rho),
107 fRadLen(0.),
108 fIntLen(0.),
109 fTemperature(0.),
110 fPressure(0.),
111 fState(kMatStateUndefined),
112 fShader(NULL),
113 fCerenkov(NULL),
114 fElement(NULL),
115 fUserExtension(0),
116 fFWExtension(0)
117{
118 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
119 fName = fName.Strip();
121 fIndex = -1;
122 fA = a;
123 fZ = z;
124 fDensity = rho;
128 SetRadLen(radlen, intlen);
129 if (!gGeoManager) {
130 gGeoManager = new TGeoManager("Geometry", "default geometry");
131 }
132 if (fZ - Int_t(fZ) > 1E-3)
133 Warning("ctor", "Material %s defined with fractional Z=%f", GetName(), fZ);
134 if (GetElement()) GetElement()->SetUsed();
136}
137
138////////////////////////////////////////////////////////////////////////////////
139/// Constructor with state, temperature and pressure.
140
142 EGeoMaterialState state, Double_t temperature, Double_t pressure)
143 :TNamed(name, ""), TAttFill(),
144 fIndex(0),
145 fA(a),
146 fZ(z),
147 fDensity(rho),
148 fRadLen(0.),
149 fIntLen(0.),
150 fTemperature(temperature),
151 fPressure(pressure),
152 fState(state),
153 fShader(NULL),
154 fCerenkov(NULL),
155 fElement(NULL),
156 fUserExtension(0),
157 fFWExtension(0)
158{
159 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
160 fName = fName.Strip();
162 fIndex = -1;
163 SetRadLen(0,0);
164 if (!gGeoManager) {
165 gGeoManager = new TGeoManager("Geometry", "default geometry");
166 }
167 if (fZ - Int_t(fZ) > 1E-3)
168 Warning("ctor", "Material %s defined with fractional Z=%f", GetName(), fZ);
169 if (GetElement()) GetElement()->SetUsed();
171}
172
173////////////////////////////////////////////////////////////////////////////////
174/// constructor
175
177 :TNamed(name, ""), TAttFill(),
178 fIndex(0),
179 fA(0.),
180 fZ(0.),
181 fDensity(rho),
182 fRadLen(0.),
183 fIntLen(0.),
184 fTemperature(0.),
185 fPressure(0.),
186 fState(kMatStateUndefined),
187 fShader(NULL),
188 fCerenkov(NULL),
189 fElement(elem),
190 fUserExtension(0),
191 fFWExtension(0)
192{
193 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
194 fName = fName.Strip();
196 fIndex = -1;
197 fA = elem->A();
198 fZ = elem->Z();
199 SetRadLen(0,0);
203 if (!gGeoManager) {
204 gGeoManager = new TGeoManager("Geometry", "default geometry");
205 }
206 if (fZ - Int_t(fZ) > 1E-3)
207 Warning("ctor", "Material %s defined with fractional Z=%f", GetName(), fZ);
208 if (GetElement()) GetElement()->SetUsed();
210}
211
212////////////////////////////////////////////////////////////////////////////////
213
215 TNamed(gm),
216 TAttFill(gm),
217 fIndex(gm.fIndex),
218 fA(gm.fA),
219 fZ(gm.fZ),
220 fDensity(gm.fDensity),
221 fRadLen(gm.fRadLen),
222 fIntLen(gm.fIntLen),
223 fTemperature(gm.fTemperature),
224 fPressure(gm.fPressure),
225 fState(gm.fState),
226 fShader(gm.fShader),
227 fCerenkov(gm.fCerenkov),
228 fElement(gm.fElement),
229 fUserExtension(gm.fUserExtension->Grab()),
230 fFWExtension(gm.fFWExtension->Grab())
231
232{
233 //copy constructor
234 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // Ensure nobody changes the units afterwards
236 TIter next(&fProperties);
237 TNamed *property;
238 while ((property = (TNamed*)next())) fProperties.Add(new TNamed(*property));
239}
240
241////////////////////////////////////////////////////////////////////////////////
242///assignment operator
243
245{
246 if(this!=&gm) {
249 fIndex=gm.fIndex;
250 fA=gm.fA;
251 fZ=gm.fZ;
253 fRadLen=gm.fRadLen;
254 fIntLen=gm.fIntLen;
257 fState=gm.fState;
258 fShader=gm.fShader;
264 TIter next(&fProperties);
265 TNamed *property;
266 while ((property = (TNamed*)next())) fProperties.Add(new TNamed(*property));
267 }
268 return *this;
269}
270
271////////////////////////////////////////////////////////////////////////////////
272/// Destructor
273
275{
278}
279
280////////////////////////////////////////////////////////////////////////////////
281/// Connect user-defined extension to the material. The material "grabs" a copy, so
282/// the original object can be released by the producer. Release the previously
283/// connected extension if any.
284///
285/// NOTE: This interface is intended for user extensions and is guaranteed not
286/// to be used by TGeo
287
289{
291 fUserExtension = 0;
292 if (ext) fUserExtension = ext->Grab();
293}
294
295//_____________________________________________________________________________
296const char *TGeoMaterial::GetPropertyRef(const char *property) const
297{
298 // Find reference for a given property
299 TNamed *prop = (TNamed*)fProperties.FindObject(property);
300 return (prop) ? prop->GetTitle() : nullptr;
301}
302
303//_____________________________________________________________________________
304TGDMLMatrix *TGeoMaterial::GetProperty(const char *property) const
305{
306 // Find reference for a given property
307 TNamed *prop = (TNamed*)fProperties.FindObject(property);
308 if ( !prop ) return nullptr;
309 return gGeoManager->GetGDMLMatrix(prop->GetTitle());
310}
311
312//_____________________________________________________________________________
314{
315 // Find reference for a given property
316 TNamed *prop = (TNamed*)fProperties.At(i);
317 if ( !prop ) return nullptr;
318 return gGeoManager->GetGDMLMatrix(prop->GetTitle());
319}
320
321//_____________________________________________________________________________
322const char *TGeoMaterial::GetConstPropertyRef(const char *property) const
323{
324 // Find reference for a given constant property
325 TNamed *prop = (TNamed*)fConstProperties.FindObject(property);
326 return (prop) ? prop->GetTitle() : nullptr;
327}
328
329//_____________________________________________________________________________
330Double_t TGeoMaterial::GetConstProperty(const char *property, Bool_t *err) const
331{
332 // Find reference for a given constant property
333 TNamed *prop = (TNamed*)fConstProperties.FindObject(property);
334 if (!prop) {
335 if (err) *err = kTRUE;
336 return 0.;
337 }
338 return gGeoManager->GetProperty(prop->GetTitle(), err);
339}
340
341//_____________________________________________________________________________
343{
344 // Find reference for a given constant property
345 TNamed *prop = (TNamed*)fConstProperties.At(i);
346 if (!prop) {
347 if (err) *err = kTRUE;
348 return 0.;
349 }
350 return gGeoManager->GetProperty(prop->GetTitle(), err);
351}
352
353//_____________________________________________________________________________
354bool TGeoMaterial::AddProperty(const char *property, const char *ref)
355{
357 if (GetPropertyRef(property)) {
358 Error("AddProperty", "Property %s already added to material %s",
359 property, GetName());
360 return false;
361 }
362 fProperties.Add(new TNamed(property, ref));
363 return true;
364}
365
366//_____________________________________________________________________________
367bool TGeoMaterial::AddConstProperty(const char *property, const char *ref)
368{
370 if (GetConstPropertyRef(property)) {
371 Error("AddConstProperty", "Constant property %s already added to material %s",
372 property, GetName());
373 return false;
374 }
375 fConstProperties.Add(new TNamed(property, ref));
376 return true;
377}
378
379////////////////////////////////////////////////////////////////////////////////
380/// Connect framework defined extension to the material. The material "grabs" a copy,
381/// so the original object can be released by the producer. Release the previously
382/// connected extension if any.
383///
384/// NOTE: This interface is intended for the use by TGeo and the users should
385/// NOT connect extensions using this method
386
388{
390 fFWExtension = 0;
391 if (ext) fFWExtension = ext->Grab();
392}
393
394////////////////////////////////////////////////////////////////////////////////
395/// Get a copy of the user extension pointer. The user must call Release() on
396/// the copy pointer once this pointer is not needed anymore (equivalent to
397/// delete() after calling new())
398
400{
401 if (fUserExtension) return fUserExtension->Grab();
402 return 0;
403}
404
405////////////////////////////////////////////////////////////////////////////////
406/// Get a copy of the framework extension pointer. The user must call Release() on
407/// the copy pointer once this pointer is not needed anymore (equivalent to
408/// delete() after calling new())
409
411{
412 if (fFWExtension) return fFWExtension->Grab();
413 return 0;
414}
415
416////////////////////////////////////////////////////////////////////////////////
417/// Provide a pointer name containing uid.
418
420{
421 static TString name;
422 name = TString::Format("pMat%d", GetUniqueID());
423 return (char*)name.Data();
424}
425
426////////////////////////////////////////////////////////////////////////////////
427/// Set radiation/absorption lengths. If the values are negative, their absolute value
428/// is taken, otherwise radlen is recomputed using G3 formula.
429
431{
432 fRadLen = TMath::Abs(radlen);
433 fIntLen = TMath::Abs(intlen);
434 // Check for vacuum
435 if (fA<0.9 || fZ<0.9) {
436 if (radlen<-1e5 || intlen<-1e-5) {
437 Error("SetRadLen","Material %s: user values taken for vacuum: radlen=%g or intlen=%g - too small", GetName(),fRadLen, fIntLen);
438 return;
439 }
440 // Ignore positive values and take big numbers
441 if (radlen>=0) fRadLen = 1.E30;
442 if (intlen>=0) fIntLen = 1.E30;
443 return;
444 }
446 // compute radlen systematically with G3 formula for a valid material
447 if ( typ == TGeoManager::kRootUnits && radlen>=0 ) {
448 //taken grom Geant3 routine GSMATE
449 constexpr Double_t alr2av = 1.39621E-03*TGeoUnit::cm2;
450 constexpr Double_t al183 = 5.20948;
454 }
455 else if ( typ == TGeoManager::kG4Units && radlen>=0 ) {
456 //taken grom Geant3 routine GSMATE
457 constexpr Double_t alr2av = 1.39621E-03*TGeant4Unit::cm2;
458 constexpr Double_t al183 = 5.20948;
462 }
463 // Compute interaction length using the same formula as in GEANT4
464 if ( typ == TGeoManager::kRootUnits && intlen>=0 ) {
465 constexpr Double_t lambda0 = 35.*TGeoUnit::g/TGeoUnit::cm2; // [g/cm^2]
466 Double_t nilinv = 0.0;
467 TGeoElement *elem = GetElement();
468 if (!elem) {
469 Fatal("SetRadLen", "Element not found for material %s", GetName());
470 return;
471 }
472 Double_t nbAtomsPerVolume = TGeoUnit::Avogadro*fDensity/elem->A();
473 nilinv += nbAtomsPerVolume*TMath::Power(elem->Neff(), 0.6666667);
474 nilinv *= TGeoUnit::amu/lambda0;
475 fIntLen = (nilinv<=0) ? TGeoShape::Big() : (TGeoUnit::cm/nilinv);
476 }
477 else if ( typ == TGeoManager::kG4Units && intlen>=0 ) {
478 constexpr Double_t lambda0 = 35.*TGeant4Unit::g/TGeant4Unit::cm2; // [g/cm^2]
479 Double_t nilinv = 0.0;
480 TGeoElement *elem = GetElement();
481 if (!elem) {
482 Fatal("SetRadLen", "Element not found for material %s", GetName());
483 return;
484 }
485 Double_t nbAtomsPerVolume = TGeant4Unit::Avogadro*fDensity/elem->A();
486 nilinv += nbAtomsPerVolume*TMath::Power(elem->Neff(), 0.6666667);
487 nilinv *= TGeant4Unit::amu/lambda0;
488 fIntLen = (nilinv<=0) ? TGeoShape::Big() : (TGeant4Unit::cm/nilinv);
489 }
490}
491
492////////////////////////////////////////////////////////////////////////////////
493/// static function
494/// Compute Coulomb correction for pair production and Brem
495/// REFERENCE : EGS MANUAL SLAC 210 - UC32 - JUNE 78
496/// FORMULA 2.7.17
497
499{
502 Double_t az2 = az*az;
503 Double_t az4 = az2 * az2;
504 Double_t fp = ( 0.0083*az4 + 0.20206 + 1./(1.+az2) ) * az2;
505 Double_t fm = ( 0.0020*az4 + 0.0369 ) * az4;
506 return fp - fm;
507}
508
509////////////////////////////////////////////////////////////////////////////////
510/// return true if the other material has the same physical properties
511
513{
514 if (other==this) return kTRUE;
515 if (other->IsMixture()) return kFALSE;
516 if (TMath::Abs(fA-other->GetA())>1E-3) return kFALSE;
517 if (TMath::Abs(fZ-other->GetZ())>1E-3) return kFALSE;
518 if (TMath::Abs(fDensity-other->GetDensity())>1E-6) return kFALSE;
519 if (GetCerenkovProperties() != other->GetCerenkovProperties()) return kFALSE;
520// if (fRadLen != other->GetRadLen()) return kFALSE;
521// if (fIntLen != other->GetIntLen()) return kFALSE;
522 return kTRUE;
523}
524
525////////////////////////////////////////////////////////////////////////////////
526/// print characteristics of this material
527
528void TGeoMaterial::Print(const Option_t * /*option*/) const
529{
530 printf("Material %s %s A=%g Z=%g rho=%g radlen=%g intlen=%g index=%i\n", GetName(), GetTitle(),
532}
533
534////////////////////////////////////////////////////////////////////////////////
535/// Save a primitive as a C++ statement(s) on output stream "out".
536
537void TGeoMaterial::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
538{
540 char *name = GetPointerName();
541 out << "// Material: " << GetName() << std::endl;
542 out << " a = " << fA << ";" << std::endl;
543 out << " z = " << fZ << ";" << std::endl;
544 out << " density = " << fDensity << ";" << std::endl;
545 out << " radl = " << fRadLen << ";" << std::endl;
546 out << " absl = " << fIntLen << ";" << std::endl;
547
548 out << " " << name << " = new TGeoMaterial(\"" << GetName() << "\", a,z,density,radl,absl);" << std::endl;
549 out << " " << name << "->SetIndex(" << GetIndex() << ");" << std::endl;
551}
552
553////////////////////////////////////////////////////////////////////////////////
554/// Get some default color related to this material.
555
557{
559 return (2+id%6);
560}
561
562////////////////////////////////////////////////////////////////////////////////
563/// Get a pointer to the element this material is made of.
564/// This second call is to avaoid warnings to not call a virtual
565/// method from the constructor
566
568{
569 if (fElement) return fElement;
571 return table->GetElement(Int_t(fZ));
572}
573
574////////////////////////////////////////////////////////////////////////////////
575/// Get a pointer to the element this material is made of.
576
578{
579 if (fElement) return fElement;
581 return table->GetElement(Int_t(fZ));
582}
583
584////////////////////////////////////////////////////////////////////////////////
585/// Single interface to get element properties.
586
588{
589 a = fA;
590 z = fZ;
591 w = 1.;
592}
593
594////////////////////////////////////////////////////////////////////////////////
595/// Retrieve material index in the list of materials
596
598{
599 if (fIndex>=0) return fIndex;
601 fIndex = matlist->IndexOf(this);
602 return fIndex;
603}
604
605////////////////////////////////////////////////////////////////////////////////
606/// Create the material representing the decay product of this material at a
607/// given time. The precision represent the minimum cumulative branching ratio for
608/// which decay products are still taken into account.
609
611{
612 TObjArray *pop = new TObjArray();
613 if (!fElement || !fElement->IsRadioNuclide()) return this;
614 FillMaterialEvolution(pop, precision);
615 Int_t ncomp = pop->GetEntriesFast();
616 if (!ncomp) return this;
617 TGeoElementRN *el;
618 Double_t *weight = new Double_t[ncomp];
619 Double_t amed = 0.;
620 Int_t i;
621 for (i=0; i<ncomp; i++) {
622 el = (TGeoElementRN *)pop->At(i);
623 weight[i] = el->Ratio()->Concentration(time) * el->A();
624 amed += weight[i];
625 }
626 Double_t rho = fDensity*amed/fA;
627 TGeoMixture *mix = 0;
628 Int_t ncomp1 = ncomp;
629 for (i=0; i<ncomp; i++) {
630 if ((weight[i]/amed)<precision) {
631 amed -= weight[i];
632 ncomp1--;
633 }
634 }
635 if (ncomp1<2) {
636 el = (TGeoElementRN *)pop->At(0);
637 delete [] weight;
638 delete pop;
639 if (ncomp1==1) return new TGeoMaterial(TString::Format("%s-evol",GetName()), el, rho);
640 return NULL;
641 }
642 mix = new TGeoMixture(TString::Format("%s-evol",GetName()), ncomp, rho);
643 for (i=0; i<ncomp; i++) {
644 weight[i] /= amed;
645 if (weight[i]<precision) continue;
646 el = (TGeoElementRN *)pop->At(i);
647 mix->AddElement(el, weight[i]);
648 }
649 delete [] weight;
650 delete pop;
651 return mix;
652}
653
654////////////////////////////////////////////////////////////////////////////////
655/// Fills a user array with all the elements deriving from the possible
656/// decay of the top element composing the mixture. Each element contained
657/// by <population> may be a radionuclide having a Bateman solution attached.
658/// The precision represent the minimum cumulative branching ratio for
659/// which decay products are still taken into account.
660/// To visualize the time evolution of each decay product one can use:
661/// ~~~ {.cpp}
662/// TGeoElement *elem = population->At(index);
663/// TGeoElementRN *elemrn = 0;
664/// if (elem->IsRadioNuclide()) elemrn = (TGeoElementRN*)elem;
665/// ~~~
666/// One can get Ni/N1(t=0) at any moment of time. Ni is the number of atoms
667/// of one of the decay products, N1(0) is the number of atoms of the top
668/// element at t=0.
669/// ~~~ {.cpp}
670/// Double_t fraction_weight = elemrn->Ratio()->Concentration(time);
671/// ~~~
672/// One can also display the time evolution of the fractional weight:
673/// ~~~ {.cpp}
674/// elemrn->Ratio()->Draw(option);
675/// ~~~
676
678{
679 if (population->GetEntriesFast()) {
680 Error("FillMaterialEvolution", "Provide an empty array !");
681 return;
682 }
684 TGeoElement *elem;
685 TGeoElementRN *elemrn;
686 TIter next(table->GetElementsRN());
687 while ((elemrn=(TGeoElementRN*)next())) elemrn->ResetRatio();
688 elem = GetElement();
689 if (!elem) {
690 Fatal("FillMaterialEvolution", "Element not found for material %s", GetName());
691 return;
692 }
693 if (!elem->IsRadioNuclide()) {
694 population->Add(elem);
695 return;
696 }
697 elemrn = (TGeoElementRN*)elem;
698 elemrn->FillPopulation(population, precision);
699}
700
701/** \class TGeoMixture
702\ingroup Geometry_classes
703
704Mixtures of elements.
705
706*/
707
709
710////////////////////////////////////////////////////////////////////////////////
711/// Default constructor
712
714{
715 fNelements = 0;
716 fZmixture = 0;
717 fAmixture = 0;
718 fWeights = 0;
719 fNatoms = 0;
721 fElements = 0;
722}
723
724////////////////////////////////////////////////////////////////////////////////
725/// constructor
726
727TGeoMixture::TGeoMixture(const char *name, Int_t /*nel*/, Double_t rho)
729{
730 fZmixture = 0;
731 fAmixture = 0;
732 fWeights = 0;
733 fNelements = 0;
734 fNatoms = 0;
736 fDensity = rho;
737 fElements = 0;
738 if (fDensity < 0) fDensity = 0.001;
739}
740
741////////////////////////////////////////////////////////////////////////////////
742/// Destructor
743
745{
746 if (fZmixture) delete[] fZmixture;
747 if (fAmixture) delete[] fAmixture;
748 if (fWeights) delete[] fWeights;
749 if (fNatoms) delete[] fNatoms;
751 if (fElements) delete fElements;
752}
753
754////////////////////////////////////////////////////////////////////////////////
755/// Compute effective A/Z and radiation length
756
758{
762 const Double_t amu = (typ==TGeoManager::kRootUnits) ? TGeoUnit::amu : TGeant4Unit::amu; // [MeV/c^2]
765 const Double_t alr2av = 1.39621E-03 * cm2;
766 const Double_t al183 = 5.20948;
767 const Double_t lambda0 = 35.*gram/cm2; // [g/cm^2]
768 Double_t radinv = 0.0;
769 Double_t nilinv = 0.0;
770 Double_t nbAtomsPerVolume;
771 fA = 0;
772 fZ = 0;
773 for (Int_t j=0;j<fNelements;j++) {
774 if (fWeights[j] <= 0) continue;
775 fA += fWeights[j]*fAmixture[j];
776 fZ += fWeights[j]*fZmixture[j];
777 nbAtomsPerVolume = na*fDensity*fWeights[j]/GetElement(j)->A();
778 nilinv += nbAtomsPerVolume*TMath::Power(GetElement(j)->Neff(), 0.6666667);
779 Double_t zc = fZmixture[j];
780 Double_t alz = TMath::Log(zc)/3.;
781 Double_t xinv = zc*(zc+TGeoMaterial::ScreenFactor(zc))*
782 (al183-alz-TGeoMaterial::Coulomb(zc))/fAmixture[j];
783 radinv += xinv*fWeights[j];
784 }
785 radinv *= alr2av*fDensity;
786 if (radinv > 0) fRadLen = cm/radinv;
787 // Compute interaction length
788 nilinv *= amu/lambda0;
789 fIntLen = (nilinv<=0) ? TGeoShape::Big() : (cm/nilinv);
790}
791
792////////////////////////////////////////////////////////////////////////////////
793/// add an element to the mixture using fraction by weight
794/// Check if the element is already defined
795
797{
799 if (z<1 || z>table->GetNelements()-1)
800 Fatal("AddElement", "Cannot add element having Z=%d to mixture %s", (Int_t)z, GetName());
801 Int_t i;
802 for (i=0; i<fNelements; i++) {
803 if (!fElements && TMath::Abs(z-fZmixture[i])<1.e-6 && TMath::Abs(a-fAmixture[i])<1.e-6) {
804 fWeights[i] += weight;
806 return;
807 }
808 }
809 if (!fNelements) {
810 fZmixture = new Double_t[1];
811 fAmixture = new Double_t[1];
812 fWeights = new Double_t[1];
813 } else {
814 Int_t nelements = fNelements+1;
815 Double_t *zmixture = new Double_t[nelements];
816 Double_t *amixture = new Double_t[nelements];
817 Double_t *weights = new Double_t[nelements];
818 for (Int_t j=0; j<fNelements; j++) {
819 zmixture[j] = fZmixture[j];
820 amixture[j] = fAmixture[j];
821 weights[j] = fWeights[j];
822 }
823 delete [] fZmixture;
824 delete [] fAmixture;
825 delete [] fWeights;
826 fZmixture = zmixture;
827 fAmixture = amixture;
828 fWeights = weights;
829 }
830
831 fNelements++;
832 i = fNelements - 1;
833 fZmixture[i] = z;
834 fAmixture[i] = a;
835 fWeights[i] = weight;
836 if (z - Int_t(z) > 1E-3)
837 Warning("DefineElement", "Mixture %s has element defined with fractional Z=%f", GetName(), z);
839 table->GetElement((Int_t)z)->SetDefined();
840
841 //compute equivalent radiation length (taken from Geant3/GSMIXT)
843}
844
845////////////////////////////////////////////////////////////////////////////////
846/// Define one component of the mixture as an existing material/mixture.
847
849{
850 TGeoElement *elnew, *elem;
851 Double_t a,z;
852 if (!mat->IsMixture()) {
853 elem = mat->GetBaseElement();
854 if (elem) {
855 AddElement(elem, weight);
856 } else {
857 a = mat->GetA();
858 z = mat->GetZ();
859 AddElement(a, z, weight);
860 }
861 return;
862 }
863 // The material is a mixture.
864 TGeoMixture *mix = (TGeoMixture*)mat;
865 Double_t wnew;
866 Int_t nelem = mix->GetNelements();
867 Bool_t elfound;
868 Int_t i,j;
869 // loop the elements of the daughter mixture
870 for (i=0; i<nelem; i++) {
871 elfound = kFALSE;
872 elnew = mix->GetElement(i);
873 if (!elnew) continue;
874 // check if we have the element already defined in the parent mixture
875 for (j=0; j<fNelements; j++) {
876 if (fWeights[j]<=0) continue;
877 elem = GetElement(j);
878 if (elem == elnew) {
879 // element found, compute new weight
880 fWeights[j] += weight * (mix->GetWmixt())[i];
881 elfound = kTRUE;
882 break;
883 }
884 }
885 if (elfound) continue;
886 // element not found, define it
887 wnew = weight * (mix->GetWmixt())[i];
888 AddElement(elnew, wnew);
889 }
890}
891
892////////////////////////////////////////////////////////////////////////////////
893/// add an element to the mixture using fraction by weight
894
896{
897 TGeoElement *elemold;
899 if (!fElements) fElements = new TObjArray(128);
900 Bool_t exist = kFALSE;
901 // If previous elements were defined by A/Z, add corresponding TGeoElements
902 for (Int_t i=0; i<fNelements; i++) {
903 elemold = (TGeoElement*)fElements->At(i);
904 if (!elemold) fElements->AddAt(elemold = table->GetElement((Int_t)fZmixture[i]), i);
905 if (elemold == elem) exist = kTRUE;
906 }
907 if (!exist) fElements->AddAtAndExpand(elem, fNelements);
908 AddElement(elem->A(), elem->Z(), weight);
909}
910
911////////////////////////////////////////////////////////////////////////////////
912/// Add a mixture element by number of atoms in the chemical formula.
913
915{
916 Int_t i,j;
917 Double_t amol;
918 TGeoElement *elemold;
920 if (!fElements) fElements = new TObjArray(128);
921 // Check if the element is already defined
922 for (i=0; i<fNelements; i++) {
923 elemold = (TGeoElement*)fElements->At(i);
924 if (!elemold) fElements->AddAt(table->GetElement((Int_t)fZmixture[i]), i);
925 else if (elemold != elem) continue;
926 if ((elem==elemold) ||
927 (TMath::Abs(elem->Z()-fZmixture[i])<1.e-6 && TMath::Abs(elem->A()-fAmixture[i])<1.e-6)) {
928 fNatoms[i] += natoms;
929 amol = 0.;
930 for (j=0; j<fNelements; j++) amol += fAmixture[j]*fNatoms[j];
931 for (j=0; j<fNelements; j++) fWeights[j] = fNatoms[j]*fAmixture[j]/amol;
933 return;
934 }
935 }
936 // New element
937 if (!fNelements) {
938 fZmixture = new Double_t[1];
939 fAmixture = new Double_t[1];
940 fWeights = new Double_t[1];
941 fNatoms = new Int_t[1];
942 } else {
943 if (!fNatoms) {
944 Fatal("AddElement", "Cannot add element by natoms in mixture %s after defining elements by weight",
945 GetName());
946 return;
947 }
948 Int_t nelements = fNelements+1;
949 Double_t *zmixture = new Double_t[nelements];
950 Double_t *amixture = new Double_t[nelements];
951 Double_t *weights = new Double_t[nelements];
952 Int_t *nnatoms = new Int_t[nelements];
953 for (j=0; j<fNelements; j++) {
954 zmixture[j] = fZmixture[j];
955 amixture[j] = fAmixture[j];
956 weights[j] = fWeights[j];
957 nnatoms[j] = fNatoms[j];
958 }
959 delete [] fZmixture;
960 delete [] fAmixture;
961 delete [] fWeights;
962 delete [] fNatoms;
963 fZmixture = zmixture;
964 fAmixture = amixture;
965 fWeights = weights;
966 fNatoms = nnatoms;
967 }
968 fNelements++;
969 Int_t iel = fNelements-1;
970 fZmixture[iel] = elem->Z();
971 fAmixture[iel] = elem->A();
972 fNatoms[iel] = natoms;
973 fElements->AddAtAndExpand(elem, iel);
974 amol = 0.;
975 for (i=0; i<fNelements; i++) {
976 if (fNatoms[i]<=0) return;
977 amol += fAmixture[i]*fNatoms[i];
978 }
979 for (i=0; i<fNelements; i++) fWeights[i] = fNatoms[i]*fAmixture[i]/amol;
980 table->GetElement(elem->Z())->SetDefined();
982}
983
984////////////////////////////////////////////////////////////////////////////////
985/// Define the mixture element at index iel by number of atoms in the chemical formula.
986
988{
990 TGeoElement *elem = table->GetElement(z);
991 if (!elem) {
992 Fatal("DefineElement", "In mixture %s, element with Z=%i not found",GetName(),z);
993 return;
994 }
995 AddElement(elem, natoms);
996}
997
998////////////////////////////////////////////////////////////////////////////////
999/// Retrieve the pointer to the element corresponding to component I.
1000
1002{
1003 if (i<0 || i>=fNelements) {
1004 Error("GetElement", "Mixture %s has only %d elements", GetName(), fNelements);
1005 return 0;
1006 }
1007 TGeoElement *elem = 0;
1008 if (fElements) elem = (TGeoElement*)fElements->At(i);
1009 if (elem) return elem;
1011 return table->GetElement(Int_t(fZmixture[i]));
1012}
1013
1014////////////////////////////////////////////////////////////////////////////////
1015/// Get specific activity (in Bq/gram) for the whole mixture (no argument) or
1016/// for a given component.
1017
1019{
1020 if (i>=0 && i<fNelements) return fWeights[i]*GetElement(i)->GetSpecificActivity();
1021 Double_t sa = 0;
1022 for (Int_t iel=0; iel<fNelements; iel++) {
1023 sa += fWeights[iel]*GetElement(iel)->GetSpecificActivity();
1024 }
1025 return sa;
1026}
1027
1028////////////////////////////////////////////////////////////////////////////////
1029/// Return true if the other material has the same physical properties
1030
1032{
1033 if (other->IsEqual(this)) return kTRUE;
1034 if (!other->IsMixture()) return kFALSE;
1035 TGeoMixture *mix = (TGeoMixture*)other;
1036 if (!mix) return kFALSE;
1037 if (fNelements != mix->GetNelements()) return kFALSE;
1038 if (TMath::Abs(fA-other->GetA())>1E-3) return kFALSE;
1039 if (TMath::Abs(fZ-other->GetZ())>1E-3) return kFALSE;
1040 if (TMath::Abs(fDensity-other->GetDensity())>1E-6) return kFALSE;
1041 if (GetCerenkovProperties() != other->GetCerenkovProperties()) return kFALSE;
1042// if (fRadLen != other->GetRadLen()) return kFALSE;
1043// if (fIntLen != other->GetIntLen()) return kFALSE;
1044 for (Int_t i=0; i<fNelements; i++) {
1045 if (TMath::Abs(fZmixture[i]-(mix->GetZmixt())[i])>1E-3) return kFALSE;
1046 if (TMath::Abs(fAmixture[i]-(mix->GetAmixt())[i])>1E-3) return kFALSE;
1047 if (TMath::Abs(fWeights[i]-(mix->GetWmixt())[i])>1E-3) return kFALSE;
1048 }
1049 return kTRUE;
1050}
1051
1052////////////////////////////////////////////////////////////////////////////////
1053/// print characteristics of this material
1054
1055void TGeoMixture::Print(const Option_t * /*option*/) const
1056{
1057 printf("Mixture %s %s Aeff=%g Zeff=%g rho=%g radlen=%g intlen=%g index=%i\n", GetName(), GetTitle(),
1059 for (Int_t i=0; i<fNelements; i++) {
1060 if (fElements && fElements->At(i)) {
1061 fElements->At(i)->Print();
1062 continue;
1063 }
1064 if (fNatoms) printf(" Element #%i : %s Z=%6.2f A=%6.2f w=%6.3f natoms=%d\n", i, GetElement(i)->GetName(),fZmixture[i],
1065 fAmixture[i], fWeights[i], fNatoms[i]);
1066 else printf(" Element #%i : %s Z=%6.2f A=%6.2f w=%6.3f\n", i, GetElement(i)->GetName(),fZmixture[i],
1067 fAmixture[i], fWeights[i]);
1068 }
1069}
1070
1071////////////////////////////////////////////////////////////////////////////////
1072/// Save a primitive as a C++ statement(s) on output stream "out".
1073
1074void TGeoMixture::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
1075{
1077 char *name = GetPointerName();
1078 out << "// Mixture: " << GetName() << std::endl;
1079 out << " nel = " << fNelements << ";" << std::endl;
1080 out << " density = " << fDensity << ";" << std::endl;
1081 out << " " << name << " = new TGeoMixture(\"" << GetName() << "\", nel,density);" << std::endl;
1082 for (Int_t i=0; i<fNelements; i++) {
1083 TGeoElement *el = GetElement(i);
1084 out << " a = " << fAmixture[i] << "; z = "<< fZmixture[i] << "; w = " << fWeights[i] << "; // " << el->GetName() << std::endl;
1085 out << " " << name << "->DefineElement(" << i << ",a,z,w);" << std::endl;
1086 }
1087 out << " " << name << "->SetIndex(" << GetIndex() << ");" << std::endl;
1089}
1090
1091////////////////////////////////////////////////////////////////////////////////
1092/// Create the mixture representing the decay product of this material at a
1093/// given time. The precision represent the minimum cumulative branching ratio for
1094/// which decay products are still taken into account.
1095
1097{
1098 TObjArray *pop = new TObjArray();
1099 FillMaterialEvolution(pop, precision);
1100 Int_t ncomp = pop->GetEntriesFast();
1101 if (!ncomp) return this;
1102 TGeoElement *elem;
1103 TGeoElementRN *el;
1104 Double_t *weight = new Double_t[ncomp];
1105 Double_t amed = 0.;
1106 Int_t i, j;
1107 for (i=0; i<ncomp; i++) {
1108 elem = (TGeoElement *)pop->At(i);
1109 if (!elem->IsRadioNuclide()) {
1110 j = fElements->IndexOf(elem);
1111 weight[i] = fWeights[j]*fAmixture[0]/fWeights[0];
1112 } else {
1113 el = (TGeoElementRN*)elem;
1114 weight[i] = el->Ratio()->Concentration(time) * el->A();
1115 }
1116 amed += weight[i];
1117 }
1118 Double_t rho = fDensity * fWeights[0] * amed/fAmixture[0];
1119 TGeoMixture *mix = 0;
1120 Int_t ncomp1 = ncomp;
1121 for (i=0; i<ncomp; i++) {
1122 if ((weight[i]/amed)<precision) {
1123 amed -= weight[i];
1124 ncomp1--;
1125 }
1126 }
1127 if (ncomp1<2) {
1128 el = (TGeoElementRN *)pop->At(0);
1129 delete [] weight;
1130 delete pop;
1131 if (ncomp1==1) return new TGeoMaterial(TString::Format("%s-evol",GetName()), el, rho);
1132 return NULL;
1133 }
1134 mix = new TGeoMixture(TString::Format("%s-evol",GetName()), ncomp, rho);
1135 for (i=0; i<ncomp; i++) {
1136 weight[i] /= amed;
1137 if (weight[i]<precision) continue;
1138 el = (TGeoElementRN *)pop->At(i);
1139 mix->AddElement(el, weight[i]);
1140 }
1141 delete [] weight;
1142 delete pop;
1143 return mix;
1144}
1145
1146////////////////////////////////////////////////////////////////////////////////
1147/// Fills a user array with all the elements deriving from the possible
1148/// decay of the top elements composing the mixture. Each element contained
1149/// by <population> may be a radionuclide having a Bateman solution attached.
1150/// The precision represent the minimum cumulative branching ratio for
1151/// which decay products are still taken into account.
1152/// To visualize the time evolution of each decay product one can use:
1153/// ~~~ {.cpp}
1154/// TGeoElement *elem = population->At(index);
1155/// TGeoElementRN *elemrn = 0;
1156/// if (elem->IsRadioNuclide()) elemrn = (TGeoElementRN*)elem;
1157/// ~~~
1158/// One can get Ni/N1(t=0) at any moment of time. Ni is the number of atoms
1159/// of one of the decay products, N1(0) is the number of atoms of the first top
1160/// element at t=0.
1161/// ~~~ {.cpp}
1162/// Double_t fraction_weight = elemrn->Ratio()->Concentration(time);
1163/// ~~~
1164/// One can also display the time evolution of the fractional weight:
1165/// ~~~ {.cpp}
1166/// elemrn->Ratio()->Draw(option);
1167/// ~~~
1168
1170{
1171 if (population->GetEntriesFast()) {
1172 Error("FillMaterialEvolution", "Provide an empty array !");
1173 return;
1174 }
1176 TGeoElement *elem;
1177 TGeoElementRN *elemrn;
1178 TIter next(table->GetElementsRN());
1179 while ((elemrn=(TGeoElementRN*)next())) elemrn->ResetRatio();
1180 Double_t factor;
1181 for (Int_t i=0; i<fNelements; i++) {
1182 elem = GetElement(i);
1183 if (!elem->IsRadioNuclide()) {
1184 population->Add(elem);
1185 continue;
1186 }
1187 elemrn = (TGeoElementRN*)elem;
1188 factor = fWeights[i]*fAmixture[0]/(fWeights[0]*fAmixture[i]);
1189 elemrn->FillPopulation(population, precision, factor);
1190 }
1191}
1192
1193////////////////////////////////////////////////////////////////////////////////
1194/// static function
1195/// Compute screening factor for pair production and Bremsstrahlung
1196/// REFERENCE : EGS MANUAL SLAC 210 - UC32 - JUNE 78
1197/// FORMULA 2.7.22
1198
1200{
1201 const Double_t al183= 5.20948 , al1440 = 7.27239;
1202 Double_t alz = TMath::Log(z)/3.;
1203 Double_t factor = (al1440 - 2*alz) / (al183 - alz - TGeoMaterial::Coulomb(z));
1204 return factor;
1205}
1206
1207////////////////////////////////////////////////////////////////////////////////
1208/// Compute Derived Quantities as in Geant4
1209
1211{
1214
1216
1218
1219 // Formula taken from G4Material.cxx L312
1220 for (Int_t i=0; i<fNelements; ++i) {
1222 }
1225}
1226
1227
1228////////////////////////////////////////////////////////////////////////////////
1229/// Compute Radiation Length based on Geant4 formula
1230
1232{
1233 // Formula taken from G4Material.cxx L556
1235 Double_t radinv = 0.0 ;
1236 for (Int_t i=0;i<fNelements;++i) {
1237 radinv += fVecNbOfAtomsPerVolume[i]*((TGeoElement*)fElements->At(i))->GetfRadTsai();
1238 }
1239 fRadLen = (radinv <= 0.0 ? DBL_MAX : cm/radinv);
1240}
1241
1242////////////////////////////////////////////////////////////////////////////////
1243/// Compute Nuclear Interaction Length based on Geant4 formula
1245{
1246 // Formula taken from G4Material.cxx L567
1251 const Double_t lambda0 = 35*g/(cm*cm);
1252 const Double_t twothird = 2.0/3.0;
1253 Double_t NILinv = 0.0;
1254 for (Int_t i=0; i<fNelements; ++i) {
1255 Int_t Z = static_cast<Int_t>(((TGeoElement*)fElements->At(i))->Z()+0.5);
1256 Double_t A = ((TGeoElement*)fElements->At(i))->Neff();
1257 if(1 == Z) {
1258 NILinv += fVecNbOfAtomsPerVolume[i]*A;
1259 } else {
1260 NILinv += fVecNbOfAtomsPerVolume[i]*TMath::Exp(twothird*TMath::Log(A));
1261 }
1262 }
1263 NILinv *= amu/lambda0;
1264 fIntLen = (NILinv <= 0.0 ? DBL_MAX : cm/NILinv);
1265}
int Int_t
Definition: CPyCppyy.h:43
#define g(i)
Definition: RSha256.hxx:105
#define e(i)
Definition: RSha256.hxx:103
const Bool_t kFALSE
Definition: RtypesCore.h:101
bool Bool_t
Definition: RtypesCore.h:63
double Double_t
Definition: RtypesCore.h:59
const Bool_t kTRUE
Definition: RtypesCore.h:100
const char Option_t
Definition: RtypesCore.h:66
#define ClassImp(name)
Definition: Rtypes.h:364
R__EXTERN TGeoManager * gGeoManager
Definition: TGeoManager.h:602
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
Int_t GetNelements() const
Definition: TGeoElement.h:417
TObjArray * GetElementsRN() const
Definition: TGeoElement.h:413
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.
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
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
TGeoElement * GetBaseElement() const
Definition: TGeoMaterial.h:113
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.
virtual TObject * GetCerenkovProperties() const
Definition: TGeoMaterial.h:118
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 * GetAmixt() const
Definition: TGeoMaterial.h:196
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 * GetZmixt() const
Definition: TGeoMaterial.h:195
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
static Double_t Big()
Definition: TGeoShape.h:88
A doubly linked list.
Definition: TList.h:44
virtual void Add(TObject *obj)
Definition: TList.h:87
virtual TObject * FindObject(const char *name) const
Find an object in this list using its name.
Definition: TList.cxx:578
virtual TObject * At(Int_t idx) const
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
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
virtual const char * GetTitle() const
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:37
Int_t IndexOf(const TObject *obj) const
Definition: TObjArray.cxx:605
Int_t GetEntriesFast() const
Definition: TObjArray.h:64
virtual void AddAtAndExpand(TObject *obj, Int_t idx)
Add object at position idx.
Definition: TObjArray.cxx:235
void Add(TObject *obj)
Definition: TObjArray.h:74
virtual void AddAt(TObject *obj, Int_t idx)
Add object at position ids.
Definition: TObjArray.cxx:254
TObject * At(Int_t idx) const
Definition: TObjArray.h:166
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 const std::string name("name")
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)
Definition: TMath.h:727
constexpr Double_t E()
Base of natural log:
Definition: TMath.h:96
Double_t Log(Double_t x)
Definition: TMath.h:760
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:735
constexpr Double_t Na()
Avogadro constant (Avogadro's Number) in .
Definition: TMath.h:291
Short_t Abs(Short_t d)
Definition: TMathBase.h:120
auto * a
Definition: textangle.C:12