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