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