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 <limits>
22#include "TMath.h"
23#include "TObjArray.h"
24#include "TGeoElement.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 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // 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 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // 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 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // 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 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // 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 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // 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 TGeoManager::SetDefaultUnits(TGeoManager::GetDefaultUnits()); // 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 == TGeoManager::kRootUnits && 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 == TGeoManager::kG4Units && 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 == TGeoManager::kRootUnits && 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 == TGeoManager::kG4Units && 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/// This second call is to avaoid warnings to not call a virtual
566/// method from the constructor
567
569{
570 if (fElement) return fElement;
572 return table->GetElement(Int_t(fZ));
573}
574
575////////////////////////////////////////////////////////////////////////////////
576/// Get a pointer to the element this material is made of.
577
579{
580 if (fElement) return fElement;
582 return table->GetElement(Int_t(fZ));
583}
584
585////////////////////////////////////////////////////////////////////////////////
586/// Single interface to get element properties.
587
589{
590 a = fA;
591 z = fZ;
592 w = 1.;
593}
594
595////////////////////////////////////////////////////////////////////////////////
596/// Retrieve material index in the list of materials
597
599{
600 if (fIndex>=0) return fIndex;
602 fIndex = matlist->IndexOf(this);
603 return fIndex;
604}
605
606////////////////////////////////////////////////////////////////////////////////
607/// Create the material representing the decay product of this material at a
608/// given time. The precision represent the minimum cumulative branching ratio for
609/// which decay products are still taken into account.
610
612{
613 TObjArray *pop = new TObjArray();
614 if (!fElement || !fElement->IsRadioNuclide()) return this;
615 FillMaterialEvolution(pop, precision);
616 Int_t ncomp = pop->GetEntriesFast();
617 if (!ncomp) return this;
618 TGeoElementRN *el;
619 Double_t *weight = new Double_t[ncomp];
620 Double_t amed = 0.;
621 Int_t i;
622 for (i=0; i<ncomp; i++) {
623 el = (TGeoElementRN *)pop->At(i);
624 weight[i] = el->Ratio()->Concentration(time) * el->A();
625 amed += weight[i];
626 }
627 Double_t rho = fDensity*amed/fA;
628 TGeoMixture *mix = 0;
629 Int_t ncomp1 = ncomp;
630 for (i=0; i<ncomp; i++) {
631 if ((weight[i]/amed)<precision) {
632 amed -= weight[i];
633 ncomp1--;
634 }
635 }
636 if (ncomp1<2) {
637 el = (TGeoElementRN *)pop->At(0);
638 delete [] weight;
639 delete pop;
640 if (ncomp1==1) return new TGeoMaterial(TString::Format("%s-evol",GetName()), el, rho);
641 return NULL;
642 }
643 mix = new TGeoMixture(TString::Format("%s-evol",GetName()), ncomp, rho);
644 for (i=0; i<ncomp; i++) {
645 weight[i] /= amed;
646 if (weight[i]<precision) continue;
647 el = (TGeoElementRN *)pop->At(i);
648 mix->AddElement(el, weight[i]);
649 }
650 delete [] weight;
651 delete pop;
652 return mix;
653}
654
655////////////////////////////////////////////////////////////////////////////////
656/// Fills a user array with all the elements deriving from the possible
657/// decay of the top element composing the mixture. Each element contained
658/// by <population> may be a radionuclide having a Bateman solution attached.
659/// The precision represent the minimum cumulative branching ratio for
660/// which decay products are still taken into account.
661/// To visualize the time evolution of each decay product one can use:
662/// ~~~ {.cpp}
663/// TGeoElement *elem = population->At(index);
664/// TGeoElementRN *elemrn = 0;
665/// if (elem->IsRadioNuclide()) elemrn = (TGeoElementRN*)elem;
666/// ~~~
667/// One can get Ni/N1(t=0) at any moment of time. Ni is the number of atoms
668/// of one of the decay products, N1(0) is the number of atoms of the top
669/// element at t=0.
670/// ~~~ {.cpp}
671/// Double_t fraction_weight = elemrn->Ratio()->Concentration(time);
672/// ~~~
673/// One can also display the time evolution of the fractional weight:
674/// ~~~ {.cpp}
675/// elemrn->Ratio()->Draw(option);
676/// ~~~
677
679{
680 if (population->GetEntriesFast()) {
681 Error("FillMaterialEvolution", "Provide an empty array !");
682 return;
683 }
685 TGeoElement *elem;
686 TGeoElementRN *elemrn;
687 TIter next(table->GetElementsRN());
688 while ((elemrn=(TGeoElementRN*)next())) elemrn->ResetRatio();
689 elem = GetElement();
690 if (!elem) {
691 Fatal("FillMaterialEvolution", "Element not found for material %s", GetName());
692 return;
693 }
694 if (!elem->IsRadioNuclide()) {
695 population->Add(elem);
696 return;
697 }
698 elemrn = (TGeoElementRN*)elem;
699 elemrn->FillPopulation(population, precision);
700}
701
702/** \class TGeoMixture
703\ingroup Geometry_classes
704
705Mixtures of elements.
706
707*/
708
710
711////////////////////////////////////////////////////////////////////////////////
712/// Default constructor
713
715{
716 fNelements = 0;
717 fZmixture = 0;
718 fAmixture = 0;
719 fWeights = 0;
720 fNatoms = 0;
722 fElements = 0;
723}
724
725////////////////////////////////////////////////////////////////////////////////
726/// constructor
727
728TGeoMixture::TGeoMixture(const char *name, Int_t /*nel*/, Double_t rho)
730{
731 fZmixture = 0;
732 fAmixture = 0;
733 fWeights = 0;
734 fNelements = 0;
735 fNatoms = 0;
737 fDensity = rho;
738 fElements = 0;
739 if (fDensity < 0) fDensity = 0.001;
740}
741
742////////////////////////////////////////////////////////////////////////////////
743/// Destructor
744
746{
747 if (fZmixture) delete[] fZmixture;
748 if (fAmixture) delete[] fAmixture;
749 if (fWeights) delete[] fWeights;
750 if (fNatoms) delete[] fNatoms;
752 if (fElements) delete fElements;
753}
754
755////////////////////////////////////////////////////////////////////////////////
756/// Compute effective A/Z and radiation length
757
759{
763 const Double_t amu = (typ==TGeoManager::kRootUnits) ? TGeoUnit::amu : TGeant4Unit::amu; // [MeV/c^2]
766 const Double_t alr2av = 1.39621E-03 * cm2;
767 const Double_t al183 = 5.20948;
768 const Double_t lambda0 = 35.*gram/cm2; // [g/cm^2]
769 Double_t radinv = 0.0;
770 Double_t nilinv = 0.0;
771 Double_t nbAtomsPerVolume;
772 fA = 0;
773 fZ = 0;
774 for (Int_t j=0;j<fNelements;j++) {
775 if (fWeights[j] <= 0) continue;
776 fA += fWeights[j]*fAmixture[j];
777 fZ += fWeights[j]*fZmixture[j];
778 nbAtomsPerVolume = na*fDensity*fWeights[j]/GetElement(j)->A();
779 nilinv += nbAtomsPerVolume*TMath::Power(GetElement(j)->Neff(), 0.6666667);
780 Double_t zc = fZmixture[j];
781 Double_t alz = TMath::Log(zc)/3.;
782 Double_t xinv = zc*(zc+TGeoMaterial::ScreenFactor(zc))*
783 (al183-alz-TGeoMaterial::Coulomb(zc))/fAmixture[j];
784 radinv += xinv*fWeights[j];
785 }
786 radinv *= alr2av*fDensity;
787 if (radinv > 0) fRadLen = cm/radinv;
788 // Compute interaction length
789 nilinv *= amu/lambda0;
790 fIntLen = (nilinv<=0) ? TGeoShape::Big() : (cm/nilinv);
791}
792
793////////////////////////////////////////////////////////////////////////////////
794/// add an element to the mixture using fraction by weight
795/// Check if the element is already defined
796
798{
800
801 // Check preconditions
802 if (weight < 0e0) {
803 Fatal("AddElement", "Cannot add element with negative weight %g to mixture %s", weight, GetName());
804 }
805 else if ( weight < std::numeric_limits<Double_t>::epsilon() ) {
806 return;
807 }
808 else if (z<1 || z>table->GetNelements()-1) {
809 Fatal("AddElement", "Cannot add element having Z=%d to mixture %s", (Int_t)z, GetName());
810 }
811 Int_t i;
812 for (i=0; i<fNelements; i++) {
813 if (!fElements && TMath::Abs(z-fZmixture[i])<1.e-6 && TMath::Abs(a-fAmixture[i])<1.e-6) {
814 fWeights[i] += weight;
816 return;
817 }
818 }
819 if (!fNelements) {
820 fZmixture = new Double_t[1];
821 fAmixture = new Double_t[1];
822 fWeights = new Double_t[1];
823 } else {
824 Int_t nelements = fNelements+1;
825 Double_t *zmixture = new Double_t[nelements];
826 Double_t *amixture = new Double_t[nelements];
827 Double_t *weights = new Double_t[nelements];
828 for (Int_t j=0; j<fNelements; j++) {
829 zmixture[j] = fZmixture[j];
830 amixture[j] = fAmixture[j];
831 weights[j] = fWeights[j];
832 }
833 delete [] fZmixture;
834 delete [] fAmixture;
835 delete [] fWeights;
836 fZmixture = zmixture;
837 fAmixture = amixture;
838 fWeights = weights;
839 }
840
841 fNelements++;
842 i = fNelements - 1;
843 fZmixture[i] = z;
844 fAmixture[i] = a;
845 fWeights[i] = weight;
846 if (z - Int_t(z) > 1E-3)
847 Warning("DefineElement", "Mixture %s has element defined with fractional Z=%f", GetName(), z);
849 table->GetElement((Int_t)z)->SetDefined();
850
851 //compute equivalent radiation length (taken from Geant3/GSMIXT)
853}
854
855////////////////////////////////////////////////////////////////////////////////
856/// Define one component of the mixture as an existing material/mixture.
857
859{
860 TGeoElement *elnew, *elem;
861 Double_t a,z;
862
863 // Check preconditions
864 if (!mat) {
865 Fatal("AddElement", "Cannot add INVALID material to mixture %s", GetName());
866 }
867 else if (weight < 0e0) {
868 Fatal("AddElement", "Cannot add material %s with negative weight %g to mixture %s",
869 mat->GetName(), weight, GetName());
870 }
871 else if ( weight < std::numeric_limits<Double_t>::epsilon() ) {
872 return;
873 }
874 if (!mat->IsMixture()) {
875 elem = mat->GetBaseElement();
876 if (elem) {
877 AddElement(elem, weight);
878 } else {
879 a = mat->GetA();
880 z = mat->GetZ();
881 AddElement(a, z, weight);
882 }
883 return;
884 }
885 // The material is a mixture.
886 TGeoMixture *mix = (TGeoMixture*)mat;
887 Double_t wnew;
888 Int_t nelem = mix->GetNelements();
889 Bool_t elfound;
890 Int_t i,j;
891 // loop the elements of the daughter mixture
892 for (i=0; i<nelem; i++) {
893 elfound = kFALSE;
894 elnew = mix->GetElement(i);
895 if (!elnew) continue;
896 // check if we have the element already defined in the parent mixture
897 for (j=0; j<fNelements; j++) {
898 if ( fWeights[j] < 0e0 ) continue;
899 elem = GetElement(j);
900 if (elem == elnew) {
901 // element found, compute new weight
902 fWeights[j] += weight * (mix->GetWmixt())[i];
903 elfound = kTRUE;
904 break;
905 }
906 }
907 if (elfound) continue;
908 // element not found, define it
909 wnew = weight * (mix->GetWmixt())[i];
910 AddElement(elnew, wnew);
911 }
912}
913
914////////////////////////////////////////////////////////////////////////////////
915/// add an element to the mixture using fraction by weight
916
918{
919 TGeoElement *elemold;
921 if (!fElements) fElements = new TObjArray(128);
922 Bool_t exist = kFALSE;
923
924 // Check preconditions
925 if (!elem) {
926 Fatal("AddElement", "Cannot add INVALID element to mixture %s", GetName());
927 }
928 else if (weight < 0e0) {
929 Fatal("AddElement", "Cannot add element %s with negative weight %g to mixture %s",
930 elem->GetName(), weight, GetName());
931 }
932 else if ( weight < std::numeric_limits<Double_t>::epsilon() ) {
933 return;
934 }
935 // If previous elements were defined by A/Z, add corresponding TGeoElements
936 for (Int_t i=0; i<fNelements; i++) {
937 elemold = (TGeoElement*)fElements->At(i);
938 if (!elemold) {
939 fElements->AddAt(elemold = table->GetElement((Int_t)fZmixture[i]), i);
940 }
941 if (elemold == elem) exist = kTRUE;
942 }
943 if (!exist) {
945 }
946 AddElement(elem->A(), elem->Z(), weight);
947}
948
949////////////////////////////////////////////////////////////////////////////////
950/// Add a mixture element by number of atoms in the chemical formula.
951
953{
954 Int_t i,j;
955 Double_t amol;
956 TGeoElement *elemold;
958 if (!fElements) fElements = new TObjArray(128);
959 // Check if the element is already defined
960 for (i=0; i<fNelements; i++) {
961 elemold = (TGeoElement*)fElements->At(i);
962 if (!elemold) fElements->AddAt(table->GetElement((Int_t)fZmixture[i]), i);
963 else if (elemold != elem) continue;
964 if ((elem==elemold) ||
965 (TMath::Abs(elem->Z()-fZmixture[i])<1.e-6 && TMath::Abs(elem->A()-fAmixture[i])<1.e-6)) {
966 fNatoms[i] += natoms;
967 amol = 0.;
968 for (j=0; j<fNelements; j++) amol += fAmixture[j]*fNatoms[j];
969 for (j=0; j<fNelements; j++) fWeights[j] = fNatoms[j]*fAmixture[j]/amol;
971 return;
972 }
973 }
974 // New element
975 if (!fNelements) {
976 fZmixture = new Double_t[1];
977 fAmixture = new Double_t[1];
978 fWeights = new Double_t[1];
979 fNatoms = new Int_t[1];
980 } else {
981 if (!fNatoms) {
982 Fatal("AddElement", "Cannot add element by natoms in mixture %s after defining elements by weight",
983 GetName());
984 return;
985 }
986 Int_t nelements = fNelements+1;
987 Double_t *zmixture = new Double_t[nelements];
988 Double_t *amixture = new Double_t[nelements];
989 Double_t *weights = new Double_t[nelements];
990 Int_t *nnatoms = new Int_t[nelements];
991 for (j=0; j<fNelements; j++) {
992 zmixture[j] = fZmixture[j];
993 amixture[j] = fAmixture[j];
994 weights[j] = fWeights[j];
995 nnatoms[j] = fNatoms[j];
996 }
997 delete [] fZmixture;
998 delete [] fAmixture;
999 delete [] fWeights;
1000 delete [] fNatoms;
1001 fZmixture = zmixture;
1002 fAmixture = amixture;
1003 fWeights = weights;
1004 fNatoms = nnatoms;
1005 }
1006 fNelements++;
1007 Int_t iel = fNelements-1;
1008 fZmixture[iel] = elem->Z();
1009 fAmixture[iel] = elem->A();
1010 fNatoms[iel] = natoms;
1011 fElements->AddAtAndExpand(elem, iel);
1012 amol = 0.;
1013 for (i=0; i<fNelements; i++) {
1014 if (fNatoms[i]<=0) return;
1015 amol += fAmixture[i]*fNatoms[i];
1016 }
1017 for (i=0; i<fNelements; i++) fWeights[i] = fNatoms[i]*fAmixture[i]/amol;
1018 table->GetElement(elem->Z())->SetDefined();
1020}
1021
1022////////////////////////////////////////////////////////////////////////////////
1023/// Define the mixture element at index iel by number of atoms in the chemical formula.
1024
1026{
1028 TGeoElement *elem = table->GetElement(z);
1029 if (!elem) {
1030 Fatal("DefineElement", "In mixture %s, element with Z=%i not found",GetName(),z);
1031 return;
1032 }
1033 AddElement(elem, natoms);
1034}
1035
1036////////////////////////////////////////////////////////////////////////////////
1037/// Retrieve the pointer to the element corresponding to component I.
1038
1040{
1041 if (i<0 || i>=fNelements) {
1042 Error("GetElement", "Mixture %s has only %d elements", GetName(), fNelements);
1043 return 0;
1044 }
1045 TGeoElement *elem = 0;
1046 if (fElements) elem = (TGeoElement*)fElements->At(i);
1047 if (elem) return elem;
1049 return table->GetElement(Int_t(fZmixture[i]));
1050}
1051
1052////////////////////////////////////////////////////////////////////////////////
1053/// Get specific activity (in Bq/gram) for the whole mixture (no argument) or
1054/// for a given component.
1055
1057{
1058 if (i>=0 && i<fNelements) return fWeights[i]*GetElement(i)->GetSpecificActivity();
1059 Double_t sa = 0;
1060 for (Int_t iel=0; iel<fNelements; iel++) {
1061 sa += fWeights[iel]*GetElement(iel)->GetSpecificActivity();
1062 }
1063 return sa;
1064}
1065
1066////////////////////////////////////////////////////////////////////////////////
1067/// Return true if the other material has the same physical properties
1068
1070{
1071 if (other->IsEqual(this)) return kTRUE;
1072 if (!other->IsMixture()) return kFALSE;
1073 TGeoMixture *mix = (TGeoMixture*)other;
1074 if (!mix) return kFALSE;
1075 if (fNelements != mix->GetNelements()) return kFALSE;
1076 if (TMath::Abs(fA-other->GetA())>1E-3) return kFALSE;
1077 if (TMath::Abs(fZ-other->GetZ())>1E-3) return kFALSE;
1078 if (TMath::Abs(fDensity-other->GetDensity())>1E-6) return kFALSE;
1079 if (GetCerenkovProperties() != other->GetCerenkovProperties()) return kFALSE;
1080// if (fRadLen != other->GetRadLen()) return kFALSE;
1081// if (fIntLen != other->GetIntLen()) return kFALSE;
1082 for (Int_t i=0; i<fNelements; i++) {
1083 if (TMath::Abs(fZmixture[i]-(mix->GetZmixt())[i])>1E-3) return kFALSE;
1084 if (TMath::Abs(fAmixture[i]-(mix->GetAmixt())[i])>1E-3) return kFALSE;
1085 if (TMath::Abs(fWeights[i]-(mix->GetWmixt())[i])>1E-3) return kFALSE;
1086 }
1087 return kTRUE;
1088}
1089
1090////////////////////////////////////////////////////////////////////////////////
1091/// print characteristics of this material
1092
1093void TGeoMixture::Print(const Option_t * /*option*/) const
1094{
1095 printf("Mixture %s %s Aeff=%g Zeff=%g rho=%g radlen=%g intlen=%g index=%i\n", GetName(), GetTitle(),
1097 for (Int_t i=0; i<fNelements; i++) {
1098 if (fElements && fElements->At(i)) {
1099 fElements->At(i)->Print();
1100 continue;
1101 }
1102 if (fNatoms) printf(" Element #%i : %s Z=%6.2f A=%6.2f w=%6.3f natoms=%d\n", i, GetElement(i)->GetName(),fZmixture[i],
1103 fAmixture[i], fWeights[i], fNatoms[i]);
1104 else printf(" Element #%i : %s Z=%6.2f A=%6.2f w=%6.3f\n", i, GetElement(i)->GetName(),fZmixture[i],
1105 fAmixture[i], fWeights[i]);
1106 }
1107}
1108
1109////////////////////////////////////////////////////////////////////////////////
1110/// Save a primitive as a C++ statement(s) on output stream "out".
1111
1112void TGeoMixture::SavePrimitive(std::ostream &out, Option_t * /*option*/ /*= ""*/)
1113{
1115 char *name = GetPointerName();
1116 out << "// Mixture: " << GetName() << std::endl;
1117 out << " nel = " << fNelements << ";" << std::endl;
1118 out << " density = " << fDensity << ";" << std::endl;
1119 out << " " << name << " = new TGeoMixture(\"" << GetName() << "\", nel,density);" << std::endl;
1120 for (Int_t i=0; i<fNelements; i++) {
1121 TGeoElement *el = GetElement(i);
1122 out << " a = " << fAmixture[i] << "; z = "<< fZmixture[i] << "; w = " << fWeights[i] << "; // " << el->GetName() << std::endl;
1123 out << " " << name << "->DefineElement(" << i << ",a,z,w);" << std::endl;
1124 }
1125 out << " " << name << "->SetIndex(" << GetIndex() << ");" << std::endl;
1127}
1128
1129////////////////////////////////////////////////////////////////////////////////
1130/// Create the mixture representing the decay product of this material at a
1131/// given time. The precision represent the minimum cumulative branching ratio for
1132/// which decay products are still taken into account.
1133
1135{
1136 TObjArray *pop = new TObjArray();
1137 FillMaterialEvolution(pop, precision);
1138 Int_t ncomp = pop->GetEntriesFast();
1139 if (!ncomp) return this;
1140 TGeoElement *elem;
1141 TGeoElementRN *el;
1142 Double_t *weight = new Double_t[ncomp];
1143 Double_t amed = 0.;
1144 Int_t i, j;
1145 for (i=0; i<ncomp; i++) {
1146 elem = (TGeoElement *)pop->At(i);
1147 if (!elem->IsRadioNuclide()) {
1148 j = fElements->IndexOf(elem);
1149 weight[i] = fWeights[j]*fAmixture[0]/fWeights[0];
1150 } else {
1151 el = (TGeoElementRN*)elem;
1152 weight[i] = el->Ratio()->Concentration(time) * el->A();
1153 }
1154 amed += weight[i];
1155 }
1156 Double_t rho = fDensity * fWeights[0] * amed/fAmixture[0];
1157 TGeoMixture *mix = 0;
1158 Int_t ncomp1 = ncomp;
1159 for (i=0; i<ncomp; i++) {
1160 if ((weight[i]/amed)<precision) {
1161 amed -= weight[i];
1162 ncomp1--;
1163 }
1164 }
1165 if (ncomp1<2) {
1166 el = (TGeoElementRN *)pop->At(0);
1167 delete [] weight;
1168 delete pop;
1169 if (ncomp1==1) return new TGeoMaterial(TString::Format("%s-evol",GetName()), el, rho);
1170 return NULL;
1171 }
1172 mix = new TGeoMixture(TString::Format("%s-evol",GetName()), ncomp, rho);
1173 for (i=0; i<ncomp; i++) {
1174 weight[i] /= amed;
1175 if (weight[i]<precision) continue;
1176 el = (TGeoElementRN *)pop->At(i);
1177 mix->AddElement(el, weight[i]);
1178 }
1179 delete [] weight;
1180 delete pop;
1181 return mix;
1182}
1183
1184////////////////////////////////////////////////////////////////////////////////
1185/// Fills a user array with all the elements deriving from the possible
1186/// decay of the top elements composing the mixture. Each element contained
1187/// by <population> may be a radionuclide having a Bateman solution attached.
1188/// The precision represent the minimum cumulative branching ratio for
1189/// which decay products are still taken into account.
1190/// To visualize the time evolution of each decay product one can use:
1191/// ~~~ {.cpp}
1192/// TGeoElement *elem = population->At(index);
1193/// TGeoElementRN *elemrn = 0;
1194/// if (elem->IsRadioNuclide()) elemrn = (TGeoElementRN*)elem;
1195/// ~~~
1196/// One can get Ni/N1(t=0) at any moment of time. Ni is the number of atoms
1197/// of one of the decay products, N1(0) is the number of atoms of the first top
1198/// element at t=0.
1199/// ~~~ {.cpp}
1200/// Double_t fraction_weight = elemrn->Ratio()->Concentration(time);
1201/// ~~~
1202/// One can also display the time evolution of the fractional weight:
1203/// ~~~ {.cpp}
1204/// elemrn->Ratio()->Draw(option);
1205/// ~~~
1206
1208{
1209 if (population->GetEntriesFast()) {
1210 Error("FillMaterialEvolution", "Provide an empty array !");
1211 return;
1212 }
1214 TGeoElement *elem;
1215 TGeoElementRN *elemrn;
1216 TIter next(table->GetElementsRN());
1217 while ((elemrn=(TGeoElementRN*)next())) elemrn->ResetRatio();
1218 Double_t factor;
1219 for (Int_t i=0; i<fNelements; i++) {
1220 elem = GetElement(i);
1221 if (!elem->IsRadioNuclide()) {
1222 population->Add(elem);
1223 continue;
1224 }
1225 elemrn = (TGeoElementRN*)elem;
1226 factor = fWeights[i]*fAmixture[0]/(fWeights[0]*fAmixture[i]);
1227 elemrn->FillPopulation(population, precision, factor);
1228 }
1229}
1230
1231////////////////////////////////////////////////////////////////////////////////
1232/// static function
1233/// Compute screening factor for pair production and Bremsstrahlung
1234/// REFERENCE : EGS MANUAL SLAC 210 - UC32 - JUNE 78
1235/// FORMULA 2.7.22
1236
1238{
1239 const Double_t al183= 5.20948 , al1440 = 7.27239;
1240 Double_t alz = TMath::Log(z)/3.;
1241 Double_t factor = (al1440 - 2*alz) / (al183 - alz - TGeoMaterial::Coulomb(z));
1242 return factor;
1243}
1244
1245////////////////////////////////////////////////////////////////////////////////
1246/// Compute Derived Quantities as in Geant4
1247
1249{
1252
1254
1256
1257 // Formula taken from G4Material.cxx L312
1258 for (Int_t i=0; i<fNelements; ++i) {
1260 }
1263}
1264
1265
1266////////////////////////////////////////////////////////////////////////////////
1267/// Compute Radiation Length based on Geant4 formula
1268
1270{
1271 // Formula taken from G4Material.cxx L556
1273 Double_t radinv = 0.0 ;
1274 for (Int_t i=0;i<fNelements;++i) {
1275 radinv += fVecNbOfAtomsPerVolume[i]*((TGeoElement*)fElements->At(i))->GetfRadTsai();
1276 }
1277 fRadLen = (radinv <= 0.0 ? DBL_MAX : cm/radinv);
1278}
1279
1280////////////////////////////////////////////////////////////////////////////////
1281/// Compute Nuclear Interaction Length based on Geant4 formula
1283{
1284 // Formula taken from G4Material.cxx L567
1289 const Double_t lambda0 = 35*g/(cm*cm);
1290 const Double_t twothird = 2.0/3.0;
1291 Double_t NILinv = 0.0;
1292 for (Int_t i=0; i<fNelements; ++i) {
1293 Int_t Z = static_cast<Int_t>(((TGeoElement*)fElements->At(i))->Z()+0.5);
1294 Double_t A = ((TGeoElement*)fElements->At(i))->Neff();
1295 if(1 == Z) {
1296 NILinv += fVecNbOfAtomsPerVolume[i]*A;
1297 } else {
1298 NILinv += fVecNbOfAtomsPerVolume[i]*TMath::Exp(twothird*TMath::Log(A));
1299 }
1300 }
1301 NILinv *= amu/lambda0;
1302 fIntLen = (NILinv <= 0.0 ? DBL_MAX : cm/NILinv);
1303}
#define g(i)
Definition: RSha256.hxx:105
#define e(i)
Definition: RSha256.hxx:103
int Int_t
Definition: RtypesCore.h:45
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
char name[80]
Definition: TGX11.cxx:110
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
TObjArray * GetElementsRN() const
Definition: TGeoElement.h:413
Int_t GetNelements() const
Definition: TGeoElement.h:417
Base class for chemical elements.
Definition: TGeoElement.h:37
virtual Double_t GetSpecificActivity() const
Definition: TGeoElement.h:84
Double_t A() const
Definition: TGeoElement.h:76
void SetDefined(Bool_t flag=kTRUE)
Definition: TGeoElement.h:90
virtual Bool_t IsRadioNuclide() const
Definition: TGeoElement.h:87
Double_t Neff() const
Returns effective number of nucleons.
Int_t Z() const
Definition: TGeoElement.h:73
void SetUsed(Bool_t flag=kTRUE)
Definition: TGeoElement.h:91
ABC for user objects attached to TGeoVolume or TGeoNode.
Definition: TGeoExtension.h:20
virtual TGeoExtension * Grab()=0
virtual void Release() const =0
The manager class for any TGeo geometry.
Definition: TGeoManager.h:45
static EDefaultUnits GetDefaultUnits()
TGeoElementTable * GetElementTable()
Returns material table. Creates it if not existing.
static void SetDefaultUnits(EDefaultUnits new_value)
Int_t AddMaterial(const TGeoMaterial *material)
Add a material to the list. Returns index of the material in list.
Double_t GetProperty(const char *name, Bool_t *error=nullptr) const
Get a user-defined property.
TGDMLMatrix * GetGDMLMatrix(const char *name) const
Get GDML matrix with a given name;.
TList * GetListOfMaterials() const
Definition: TGeoManager.h:491
Base class describing materials.
Definition: TGeoMaterial.h:36
virtual ~TGeoMaterial()
Destructor.
Double_t GetConstProperty(const char *property, Bool_t *error=nullptr) const
void SetUserExtension(TGeoExtension *ext)
Connect user-defined extension to the material.
char * GetPointerName() const
Provide a pointer name containing uid.
virtual TObject * GetCerenkovProperties() const
Definition: TGeoMaterial.h:118
EGeoMaterialState fState
Definition: TGeoMaterial.h:58
static Double_t ScreenFactor(Double_t z)
static function Compute screening factor for pair production and Bremsstrahlung REFERENCE : EGS MANUA...
const char * GetConstPropertyRef(const char *property) const
virtual Bool_t IsMixture() const
Definition: TGeoMaterial.h:130
Double_t fPressure
Definition: TGeoMaterial.h:57
void SetFWExtension(TGeoExtension *ext)
Connect framework defined extension to the material.
bool AddConstProperty(const char *property, const char *ref)
Double_t fTemperature
Definition: TGeoMaterial.h:56
virtual void GetElementProp(Double_t &a, Double_t &z, Double_t &w, Int_t i=0)
Single interface to get element properties.
virtual void Print(const Option_t *option="") const
print characteristics of this material
virtual TGeoMaterial * DecayMaterial(Double_t time, Double_t precision=0.001)
Create the material representing the decay product of this material at a given time.
TList fProperties
Definition: TGeoMaterial.h:62
bool AddProperty(const char *property, const char *ref)
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
Double_t fZ
Definition: TGeoMaterial.h:52
void SetRadLen(Double_t radlen, Double_t intlen=0.)
Set radiation/absorption lengths.
TGeoElement * GetElement() const
Get a pointer to the element this material is made of.
virtual Bool_t IsEq(const TGeoMaterial *other) const
return true if the other material has the same physical properties
TGeoElement * GetBaseElement() const
Definition: TGeoMaterial.h:113
const char * GetPropertyRef(const char *property) const
Double_t fA
Definition: TGeoMaterial.h:51
TGDMLMatrix * GetProperty(const char *name) const
TObject * fCerenkov
Definition: TGeoMaterial.h:60
Double_t fDensity
Definition: TGeoMaterial.h:53
void SetUsed(Bool_t flag=kTRUE)
Definition: TGeoMaterial.h:139
virtual void FillMaterialEvolution(TObjArray *population, Double_t precision=0.001)
Fills a user array with all the elements deriving from the possible decay of the top element composin...
Int_t GetIndex()
Retrieve material index in the list of materials.
Double_t fIntLen
Definition: TGeoMaterial.h:55
TGeoExtension * fUserExtension
Definition: TGeoMaterial.h:64
TList fConstProperties
Definition: TGeoMaterial.h:63
TGeoElement * fElement
Definition: TGeoMaterial.h:61
TGeoMaterial & operator=(const TGeoMaterial &)
assignment operator
TObject * fShader
Definition: TGeoMaterial.h:59
TGeoExtension * GrabUserExtension() const
Get a copy of the user extension pointer.
Double_t fRadLen
Definition: TGeoMaterial.h:54
TGeoExtension * GrabFWExtension() const
Get a copy of the framework extension pointer.
virtual Int_t GetDefaultColor() const
Get some default color related to this material.
TGeoMaterial()
Default constructor.
static Double_t Coulomb(Double_t z)
static function Compute Coulomb correction for pair production and Brem REFERENCE : EGS MANUAL SLAC 2...
virtual Double_t GetA() const
Definition: TGeoMaterial.h:105
TGeoExtension * fFWExtension
Transient user-defined extension to materials.
Definition: TGeoMaterial.h:65
virtual Double_t GetDensity() const
Definition: TGeoMaterial.h:108
virtual Double_t GetZ() const
Definition: TGeoMaterial.h:106
Mixtures of elements.
Definition: TGeoMaterial.h:157
TObjArray * fElements
Definition: TGeoMaterial.h:166
virtual TGeoElement * GetElement(Int_t i=0) const
Retrieve the pointer to the element corresponding to component I.
void ComputeNuclearInterLength()
Compute Nuclear Interaction Length based on Geant4 formula.
Double_t * GetZmixt() const
Definition: TGeoMaterial.h:195
void AddElement(Double_t a, Double_t z, Double_t weight)
add an element to the mixture using fraction by weight Check if the element is already defined
Double_t * fZmixture
Definition: TGeoMaterial.h:161
virtual Double_t GetSpecificActivity(Int_t i=-1) const
Get specific activity (in Bq/gram) for the whole mixture (no argument) or for a given component.
TGeoMixture()
Default constructor.
virtual ~TGeoMixture()
Destructor.
void ComputeDerivedQuantities()
Compute Derived Quantities as in Geant4.
virtual void Print(const Option_t *option="") const
print characteristics of this material
void ComputeRadiationLength()
Compute Radiation Length based on Geant4 formula.
Double_t * fVecNbOfAtomsPerVolume
Definition: TGeoMaterial.h:165
virtual void FillMaterialEvolution(TObjArray *population, Double_t precision=0.001)
Fills a user array with all the elements deriving from the possible decay of the top elements composi...
Double_t * fAmixture
Definition: TGeoMaterial.h:162
void AverageProperties()
Compute effective A/Z and radiation length.
Int_t fNelements
Definition: TGeoMaterial.h:160
Double_t * GetWmixt() const
Definition: TGeoMaterial.h:197
virtual TGeoMaterial * DecayMaterial(Double_t time, Double_t precision=0.001)
Create the mixture representing the decay product of this material at a given time.
virtual void SavePrimitive(std::ostream &out, Option_t *option="")
Save a primitive as a C++ statement(s) on output stream "out".
virtual Int_t GetNelements() const
Definition: TGeoMaterial.h:194
Int_t * fNatoms
Definition: TGeoMaterial.h:164
virtual Bool_t IsEq(const TGeoMaterial *other) const
Return true if the other material has the same physical properties.
Double_t * fWeights
Definition: TGeoMaterial.h:163
void DefineElement(Int_t iel, Double_t a, Double_t z, Double_t weight)
Definition: TGeoMaterial.h:215
Double_t * GetAmixt() const
Definition: TGeoMaterial.h:196
static Double_t Big()
Definition: TGeoShape.h:88
A doubly linked list.
Definition: TList.h: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
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: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 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
REAL epsilon
Definition: triangle.c:618