Geant4 11.2.2
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4MaterialPropertiesTable Class Reference

#include <G4MaterialPropertiesTable.hh>

+ Inheritance diagram for G4MaterialPropertiesTable:

Public Member Functions

 G4MaterialPropertiesTable ()
 
virtual ~G4MaterialPropertiesTable ()
 
void AddConstProperty (const G4String &key, G4double propertyValue, G4bool createNewKey=false)
 
void AddConstProperty (const char *key, G4double propertyValue, G4bool createNewKey=false)
 
G4MaterialPropertyVectorAddProperty (const G4String &key, const std::vector< G4double > &photonEnergies, const std::vector< G4double > &propertyValues, G4bool createNewKey=false, G4bool spline=false)
 
G4MaterialPropertyVectorAddProperty (const char *key, G4double *photonEnergies, G4double *propertyValues, G4int numEntries, G4bool createNewKey=false, G4bool spline=false)
 
void AddProperty (const G4String &key, G4MaterialPropertyVector *opv, G4bool createNewKey=false)
 
void AddProperty (const char *key, G4MaterialPropertyVector *opv, G4bool createNewKey=false)
 
void AddProperty (const G4String &key, const G4String &mat)
 
void RemoveConstProperty (const G4String &key)
 
void RemoveConstProperty (const char *key)
 
void RemoveProperty (const G4String &key)
 
void RemoveProperty (const char *key)
 
G4double GetConstProperty (const G4String &key) const
 
G4double GetConstProperty (const char *key) const
 
G4double GetConstProperty (const G4int index) const
 
G4bool ConstPropertyExists (const G4String &key) const
 
G4bool ConstPropertyExists (const char *key) const
 
G4bool ConstPropertyExists (const G4int index) const
 
G4MaterialPropertyVectorGetProperty (const char *key) const
 
G4MaterialPropertyVectorGetProperty (const G4String &key) const
 
G4MaterialPropertyVectorGetProperty (const G4int index) const
 
void AddEntry (const G4String &key, G4double aPhotonEnergy, G4double aPropertyValue)
 
void AddEntry (const char *key, G4double aPhotonEnergy, G4double aPropertyValue)
 
G4int GetConstPropertyIndex (const G4String &key) const
 
G4int GetPropertyIndex (const G4String &key) const
 
void DumpTable () const
 
const std::vector< G4String > & GetMaterialPropertyNames () const
 
const std::vector< G4String > & GetMaterialConstPropertyNames () const
 
const std::vector< G4MaterialPropertyVector * > & GetProperties () const
 
const std::vector< std::pair< G4double, G4bool > > & GetConstProperties () const
 

Detailed Description

Definition at line 59 of file G4MaterialPropertiesTable.hh.

Constructor & Destructor Documentation

◆ G4MaterialPropertiesTable()

G4MaterialPropertiesTable::G4MaterialPropertiesTable ( )

Definition at line 62 of file G4MaterialPropertiesTable.cc.

63{
64 // elements of these 2 vectors must be in same order as
65 // the corresponding enums in G4MaterialPropertiesIndex.hh
66 fMatPropNames.assign(kNumberOfPropertyIndex, "");
67 fMatPropNames[kRINDEX] = "RINDEX";
68 fMatPropNames[kREFLECTIVITY] = "REFLECTIVITY";
69 fMatPropNames[kREALRINDEX] = "REALRINDEX";
70 fMatPropNames[kIMAGINARYRINDEX] = "IMAGINARYRINDEX";
71 fMatPropNames[kEFFICIENCY] = "EFFICIENCY";
72 fMatPropNames[kTRANSMITTANCE] = "TRANSMITTANCE";
73 fMatPropNames[kSPECULARLOBECONSTANT] = "SPECULARLOBECONSTANT";
74 fMatPropNames[kSPECULARSPIKECONSTANT] = "SPECULARSPIKECONSTANT";
75 fMatPropNames[kBACKSCATTERCONSTANT] = "BACKSCATTERCONSTANT";
76 fMatPropNames[kGROUPVEL] = "GROUPVEL";
77 fMatPropNames[kMIEHG] = "MIEHG";
78 fMatPropNames[kRAYLEIGH] = "RAYLEIGH";
79 fMatPropNames[kWLSCOMPONENT] = "WLSCOMPONENT";
80 fMatPropNames[kWLSABSLENGTH] = "WLSABSLENGTH";
81 fMatPropNames[kWLSCOMPONENT2] = "WLSCOMPONENT2";
82 fMatPropNames[kWLSABSLENGTH2] = "WLSABSLENGTH2";
83 fMatPropNames[kABSLENGTH] = "ABSLENGTH";
84 fMatPropNames[kPROTONSCINTILLATIONYIELD] = "PROTONSCINTILLATIONYIELD";
85 fMatPropNames[kDEUTERONSCINTILLATIONYIELD] ="DEUTERONSCINTILLATIONYIELD";
86 fMatPropNames[kTRITONSCINTILLATIONYIELD] = "TRITONSCINTILLATIONYIELD";
87 fMatPropNames[kALPHASCINTILLATIONYIELD] = "ALPHASCINTILLATIONYIELD";
88 fMatPropNames[kIONSCINTILLATIONYIELD] = "IONSCINTILLATIONYIELD";
89 fMatPropNames[kELECTRONSCINTILLATIONYIELD] ="ELECTRONSCINTILLATIONYIELD";
90 fMatPropNames[kSCINTILLATIONCOMPONENT1] = "SCINTILLATIONCOMPONENT1";
91 fMatPropNames[kSCINTILLATIONCOMPONENT2] = "SCINTILLATIONCOMPONENT2";
92 fMatPropNames[kSCINTILLATIONCOMPONENT3] = "SCINTILLATIONCOMPONENT3";
93 fMatPropNames[kCOATEDRINDEX] = "COATEDRINDEX";
94
95 fMP.assign(kNumberOfPropertyIndex, nullptr);
96
97 fMatConstPropNames.assign(kNumberOfConstPropertyIndex, "");
98 fMatConstPropNames[kSURFACEROUGHNESS] = "SURFACEROUGHNESS";
99 fMatConstPropNames[kISOTHERMAL_COMPRESSIBILITY] = "ISOTHERMAL_COMPRESSIBILITY";
100 fMatConstPropNames[kRS_SCALE_FACTOR] = "RS_SCALE_FACTOR";
101 fMatConstPropNames[kWLSMEANNUMBERPHOTONS] = "WLSMEANNUMBERPHOTONS";
102 fMatConstPropNames[kWLSTIMECONSTANT] = "WLSTIMECONSTANT";
103 fMatConstPropNames[kWLSMEANNUMBERPHOTONS2] = "WLSMEANNUMBERPHOTONS2";
104 fMatConstPropNames[kWLSTIMECONSTANT2] = "WLSTIMECONSTANT2";
105 fMatConstPropNames[kMIEHG_FORWARD] = "MIEHG_FORWARD";
106 fMatConstPropNames[kMIEHG_BACKWARD] = "MIEHG_BACKWARD";
107 fMatConstPropNames[kMIEHG_FORWARD_RATIO] = "MIEHG_FORWARD_RATIO";
108 fMatConstPropNames[kSCINTILLATIONYIELD] = "SCINTILLATIONYIELD";
109 fMatConstPropNames[kRESOLUTIONSCALE] = "RESOLUTIONSCALE";
110 fMatConstPropNames[kFERMIPOT] = "FERMIPOT";
111 fMatConstPropNames[kDIFFUSION] = "DIFFUSION";
112 fMatConstPropNames[kSPINFLIP] = "SPINFLIP";
113 fMatConstPropNames[kLOSS] = "LOSS";
114 fMatConstPropNames[kLOSSCS] = "LOSSCS";
115 fMatConstPropNames[kABSCS] = "ABSCS";
116 fMatConstPropNames[kSCATCS] = "SCATCS";
117 fMatConstPropNames[kMR_NBTHETA] = "MR_NBTHETA";
118 fMatConstPropNames[kMR_NBE] = "MR_NBE";
119 fMatConstPropNames[kMR_RRMS] = "MR_RRMS";
120 fMatConstPropNames[kMR_CORRLEN] = "MR_CORRLEN";
121 fMatConstPropNames[kMR_THETAMIN] = "MR_THETAMIN";
122 fMatConstPropNames[kMR_THETAMAX] = "MR_THETAMAX";
123 fMatConstPropNames[kMR_EMIN] = "MR_EMIN";
124 fMatConstPropNames[kMR_EMAX] = "MR_EMAX";
125 fMatConstPropNames[kMR_ANGNOTHETA] = "MR_ANGNOTHETA";
126 fMatConstPropNames[kMR_ANGNOPHI] = "MR_ANGNOPHI";
127 fMatConstPropNames[kMR_ANGCUT] = "MR_ANGCUT";
128 fMatConstPropNames[kSCINTILLATIONTIMECONSTANT1] = "SCINTILLATIONTIMECONSTANT1";
129 fMatConstPropNames[kSCINTILLATIONTIMECONSTANT2] = "SCINTILLATIONTIMECONSTANT2";
130 fMatConstPropNames[kSCINTILLATIONTIMECONSTANT3] = "SCINTILLATIONTIMECONSTANT3";
131 fMatConstPropNames[kSCINTILLATIONRISETIME1] = "SCINTILLATIONRISETIME1";
132 fMatConstPropNames[kSCINTILLATIONRISETIME2] = "SCINTILLATIONRISETIME2";
133 fMatConstPropNames[kSCINTILLATIONRISETIME3] = "SCINTILLATIONRISETIME3";
134 fMatConstPropNames[kSCINTILLATIONYIELD1] = "SCINTILLATIONYIELD1";
135 fMatConstPropNames[kSCINTILLATIONYIELD2] = "SCINTILLATIONYIELD2";
136 fMatConstPropNames[kSCINTILLATIONYIELD3] = "SCINTILLATIONYIELD3";
137 fMatConstPropNames[kPROTONSCINTILLATIONYIELD1] = "PROTONSCINTILLATIONYIELD1";
138 fMatConstPropNames[kPROTONSCINTILLATIONYIELD2] = "PROTONSCINTILLATIONYIELD2";
139 fMatConstPropNames[kPROTONSCINTILLATIONYIELD3] = "PROTONSCINTILLATIONYIELD3";
140 fMatConstPropNames[kDEUTERONSCINTILLATIONYIELD1] = "DEUTERONSCINTILLATIONYIELD1";
141 fMatConstPropNames[kDEUTERONSCINTILLATIONYIELD2] = "DEUTERONSCINTILLATIONYIELD2";
142 fMatConstPropNames[kDEUTERONSCINTILLATIONYIELD3] = "DEUTERONSCINTILLATIONYIELD3";
143 fMatConstPropNames[kTRITONSCINTILLATIONYIELD1] = "TRITONSCINTILLATIONYIELD1";
144 fMatConstPropNames[kTRITONSCINTILLATIONYIELD2] = "TRITONSCINTILLATIONYIELD2";
145 fMatConstPropNames[kTRITONSCINTILLATIONYIELD3] = "TRITONSCINTILLATIONYIELD3";
146 fMatConstPropNames[kALPHASCINTILLATIONYIELD1] = "ALPHASCINTILLATIONYIELD1";
147 fMatConstPropNames[kALPHASCINTILLATIONYIELD2] = "ALPHASCINTILLATIONYIELD2";
148 fMatConstPropNames[kALPHASCINTILLATIONYIELD3] = "ALPHASCINTILLATIONYIELD3";
149 fMatConstPropNames[kIONSCINTILLATIONYIELD1] = "IONSCINTILLATIONYIELD1";
150 fMatConstPropNames[kIONSCINTILLATIONYIELD2] = "IONSCINTILLATIONYIELD2";
151 fMatConstPropNames[kIONSCINTILLATIONYIELD3] = "IONSCINTILLATIONYIELD3";
152 fMatConstPropNames[kELECTRONSCINTILLATIONYIELD1] = "ELECTRONSCINTILLATIONYIELD1";
153 fMatConstPropNames[kELECTRONSCINTILLATIONYIELD2] = "ELECTRONSCINTILLATIONYIELD2";
154 fMatConstPropNames[kELECTRONSCINTILLATIONYIELD3] = "ELECTRONSCINTILLATIONYIELD3";
155 fMatConstPropNames[kCOATEDTHICKNESS] = "COATEDTHICKNESS";
156 fMatConstPropNames[kCOATEDFRUSTRATEDTRANSMISSION] = "COATEDFRUSTRATEDTRANSMISSION";
157 fMatConstPropNames[kPROTONSCINTILLATIONTIMECONSTANT1] = "PROTONSCINTILLATIONTIMECONSTANT1";
158 fMatConstPropNames[kPROTONSCINTILLATIONTIMECONSTANT2] = "PROTONSCINTILLATIONTIMECONSTANT2";
159 fMatConstPropNames[kPROTONSCINTILLATIONTIMECONSTANT3] = "PROTONSCINTILLATIONTIMECONSTANT3";
160 fMatConstPropNames[kDEUTERONSCINTILLATIONTIMECONSTANT1] = "DEUTERONSCINTILLATIONTIMECONSTANT1";
161 fMatConstPropNames[kDEUTERONSCINTILLATIONTIMECONSTANT2] = "DEUTERONSCINTILLATIONTIMECONSTANT2";
162 fMatConstPropNames[kDEUTERONSCINTILLATIONTIMECONSTANT3] = "DEUTERONSCINTILLATIONTIMECONSTANT3";
163 fMatConstPropNames[kTRITONSCINTILLATIONTIMECONSTANT1] = "TRITONSCINTILLATIONTIMECONSTANT1";
164 fMatConstPropNames[kTRITONSCINTILLATIONTIMECONSTANT2] = "TRITONSCINTILLATIONTIMECONSTANT2";
165 fMatConstPropNames[kTRITONSCINTILLATIONTIMECONSTANT3] = "TRITONSCINTILLATIONTIMECONSTANT3";
166 fMatConstPropNames[kALPHASCINTILLATIONTIMECONSTANT1] = "ALPHASCINTILLATIONTIMECONSTANT1";
167 fMatConstPropNames[kALPHASCINTILLATIONTIMECONSTANT2] = "ALPHASCINTILLATIONTIMECONSTANT2";
168 fMatConstPropNames[kALPHASCINTILLATIONTIMECONSTANT3] = "ALPHASCINTILLATIONTIMECONSTANT3";
169 fMatConstPropNames[kIONSCINTILLATIONTIMECONSTANT1] = "IONSCINTILLATIONTIMECONSTANT1";
170 fMatConstPropNames[kIONSCINTILLATIONTIMECONSTANT2] = "IONSCINTILLATIONTIMECONSTANT2";
171 fMatConstPropNames[kIONSCINTILLATIONTIMECONSTANT3] = "IONSCINTILLATIONTIMECONSTANT3";
172 fMatConstPropNames[kELECTRONSCINTILLATIONTIMECONSTANT1] = "ELECTRONSCINTILLATIONTIMECONSTANT1";
173 fMatConstPropNames[kELECTRONSCINTILLATIONTIMECONSTANT2] = "ELECTRONSCINTILLATIONTIMECONSTANT2";
174 fMatConstPropNames[kELECTRONSCINTILLATIONTIMECONSTANT3] = "ELECTRONSCINTILLATIONTIMECONSTANT3";
175 fMCP.assign(kNumberOfConstPropertyIndex, {0., false});
176}
@ kSCINTILLATIONCOMPONENT1
@ kSCINTILLATIONCOMPONENT2
@ kSCINTILLATIONCOMPONENT3
@ kELECTRONSCINTILLATIONYIELD
@ kALPHASCINTILLATIONYIELD
@ kPROTONSCINTILLATIONYIELD
@ kDEUTERONSCINTILLATIONYIELD
@ kTRITONSCINTILLATIONYIELD
@ kSCINTILLATIONTIMECONSTANT1
@ kTRITONSCINTILLATIONYIELD1
@ kDEUTERONSCINTILLATIONYIELD3
@ kPROTONSCINTILLATIONTIMECONSTANT2
@ kALPHASCINTILLATIONTIMECONSTANT1
@ kELECTRONSCINTILLATIONTIMECONSTANT2
@ kDEUTERONSCINTILLATIONYIELD2
@ kELECTRONSCINTILLATIONTIMECONSTANT3
@ kDEUTERONSCINTILLATIONTIMECONSTANT1
@ kTRITONSCINTILLATIONTIMECONSTANT2
@ kTRITONSCINTILLATIONYIELD2
@ kNumberOfConstPropertyIndex
@ kALPHASCINTILLATIONYIELD2
@ kELECTRONSCINTILLATIONYIELD3
@ kTRITONSCINTILLATIONTIMECONSTANT1
@ kALPHASCINTILLATIONYIELD1
@ kALPHASCINTILLATIONTIMECONSTANT2
@ kPROTONSCINTILLATIONTIMECONSTANT3
@ kDEUTERONSCINTILLATIONTIMECONSTANT3
@ kELECTRONSCINTILLATIONYIELD2
@ kPROTONSCINTILLATIONYIELD2
@ kDEUTERONSCINTILLATIONYIELD1
@ kTRITONSCINTILLATIONTIMECONSTANT3
@ kISOTHERMAL_COMPRESSIBILITY
@ kTRITONSCINTILLATIONYIELD3
@ kSCINTILLATIONTIMECONSTANT3
@ kIONSCINTILLATIONTIMECONSTANT1
@ kIONSCINTILLATIONTIMECONSTANT3
@ kPROTONSCINTILLATIONYIELD3
@ kIONSCINTILLATIONTIMECONSTANT2
@ kALPHASCINTILLATIONTIMECONSTANT3
@ kELECTRONSCINTILLATIONYIELD1
@ kALPHASCINTILLATIONYIELD3
@ kSCINTILLATIONTIMECONSTANT2
@ kPROTONSCINTILLATIONYIELD1
@ kCOATEDFRUSTRATEDTRANSMISSION
@ kDEUTERONSCINTILLATIONTIMECONSTANT2
@ kPROTONSCINTILLATIONTIMECONSTANT1
@ kELECTRONSCINTILLATIONTIMECONSTANT1

◆ ~G4MaterialPropertiesTable()

G4MaterialPropertiesTable::~G4MaterialPropertiesTable ( )
virtual

Definition at line 178 of file G4MaterialPropertiesTable.cc.

179{
180 for (auto prop : fMP) {
181 delete (prop);
182 }
183}

Member Function Documentation

◆ AddConstProperty() [1/2]

void G4MaterialPropertiesTable::AddConstProperty ( const char * key,
G4double propertyValue,
G4bool createNewKey = false )

Definition at line 459 of file G4MaterialPropertiesTable.cc.

461{
462 // Provides a way of adding a constant property to the Material Properties
463 // Table given a key
464 AddConstProperty(G4String(key), propertyValue, createNewKey);
465}
void AddConstProperty(const G4String &key, G4double propertyValue, G4bool createNewKey=false)

◆ AddConstProperty() [2/2]

void G4MaterialPropertiesTable::AddConstProperty ( const G4String & key,
G4double propertyValue,
G4bool createNewKey = false )

Definition at line 434 of file G4MaterialPropertiesTable.cc.

436{
437 // Provides a way of adding a constant property to the Material Properties
438 // Table given a key
439 if (std::find(fMatConstPropNames.cbegin(), fMatConstPropNames.cend(), key) ==
440 fMatConstPropNames.cend())
441 {
442 if (createNewKey) {
443 fMatConstPropNames.push_back(key);
444 fMCP.emplace_back(0., true);
445 }
446 else {
448 ed << "Attempting to create a new material constant property key " << key
449 << " without setting\n"
450 << "createNewKey parameter of AddProperty to true.";
451 G4Exception("G4MaterialPropertiesTable::AddProperty()", "mat207", FatalException, ed);
452 }
453 }
454 G4int index = GetConstPropertyIndex(key);
455
456 fMCP[index] = std::pair<G4double, G4bool>{propertyValue, true};
457}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
std::ostringstream G4ExceptionDescription
int G4int
Definition G4Types.hh:85
G4int GetConstPropertyIndex(const G4String &key) const

Referenced by AddConstProperty(), G4GDMLReadMaterials::PropertyRead(), G4GDMLReadSolids::PropertyRead(), and G4UCNMaterialPropertiesTable::SetMicroRoughnessParameters().

◆ AddEntry() [1/2]

void G4MaterialPropertiesTable::AddEntry ( const char * key,
G4double aPhotonEnergy,
G4double aPropertyValue )

Definition at line 525 of file G4MaterialPropertiesTable.cc.

527{
528 AddEntry(G4String(key), aPhotonEnergy, aPropertyValue);
529}
void AddEntry(const G4String &key, G4double aPhotonEnergy, G4double aPropertyValue)

◆ AddEntry() [2/2]

void G4MaterialPropertiesTable::AddEntry ( const G4String & key,
G4double aPhotonEnergy,
G4double aPropertyValue )

Definition at line 489 of file G4MaterialPropertiesTable.cc.

491{
492 // Allows to add an entry pair directly to the Material Property Vector
493 // given a key.
494 if (std::find(fMatPropNames.cbegin(), fMatPropNames.cend(), key) == fMatPropNames.cend()) {
496 ed << "Material Property Vector " << key << " not found.";
497 G4Exception("G4MaterialPropertiesTable::AddEntry()", "mat214", FatalException, ed);
498 }
499 G4int index = GetPropertyIndex(key);
500
501 G4MaterialPropertyVector* targetVector = fMP[index];
502 if (targetVector != nullptr) {
503 // do not allow duplicate energies
504 for (std::size_t i = 0; i < targetVector->GetVectorLength(); ++i) {
505 if (aPhotonEnergy == targetVector->Energy(i)) {
507 ed << "Energy values in material property vector must be unique. "
508 << "Key: " << key;
509 G4Exception("G4MaterialPropertiesTable::AddEntry()", "mat217", FatalException, ed);
510 }
511 }
512
513 targetVector->InsertValues(aPhotonEnergy, aPropertyValue);
514 }
515 else {
517 ed << "Material Property Vector " << key << " not found.";
518 G4Exception("G4MaterialPropertiesTable::AddEntry()", "mat208", FatalException, ed);
519 }
520 if (key == "RINDEX") {
521 CalculateGROUPVEL();
522 }
523}
G4int GetPropertyIndex(const G4String &key) const
void InsertValues(const G4double energy, const G4double value)
G4double Energy(const std::size_t index) const
std::size_t GetVectorLength() const

Referenced by AddEntry().

◆ AddProperty() [1/5]

G4MaterialPropertyVector * G4MaterialPropertiesTable::AddProperty ( const char * key,
G4double * photonEnergies,
G4double * propertyValues,
G4int numEntries,
G4bool createNewKey = false,
G4bool spline = false )

Definition at line 359 of file G4MaterialPropertiesTable.cc.

362{
363 // Provides a way of adding a property to the Material Properties
364 // Table given a pair of arrays and a key
365 G4String k(key);
366
367 std::vector<G4double> energies(photonEnergies, photonEnergies + numEntries);
368 std::vector<G4double> values(propertyValues, propertyValues + numEntries);
369 return AddProperty(k, energies, values, createNewKey, spline);
370}
G4MaterialPropertyVector * AddProperty(const G4String &key, const std::vector< G4double > &photonEnergies, const std::vector< G4double > &propertyValues, G4bool createNewKey=false, G4bool spline=false)

◆ AddProperty() [2/5]

void G4MaterialPropertiesTable::AddProperty ( const char * key,
G4MaterialPropertyVector * opv,
G4bool createNewKey = false )

Definition at line 421 of file G4MaterialPropertiesTable.cc.

423{
424 AddProperty(G4String(key), mpv, createNewKey);
425}

◆ AddProperty() [3/5]

void G4MaterialPropertiesTable::AddProperty ( const G4String & key,
const G4String & mat )

Definition at line 427 of file G4MaterialPropertiesTable.cc.

428{
429 // load a material property vector defined in Geant4 source
431 AddProperty(key, v);
432}
G4MaterialPropertyVector * GetProperty(const G4String &key, const G4String &mat)

◆ AddProperty() [4/5]

G4MaterialPropertyVector * G4MaterialPropertiesTable::AddProperty ( const G4String & key,
const std::vector< G4double > & photonEnergies,
const std::vector< G4double > & propertyValues,
G4bool createNewKey = false,
G4bool spline = false )

Definition at line 299 of file G4MaterialPropertiesTable.cc.

302{
303 if (photonEnergies.size() != propertyValues.size()) {
305 ed << "AddProperty error. Number of property values must be equal to the number of\n"
306 << "energy values. Property name: " << key;
307 G4Exception("G4MaterialPropertiesTable::AddProperty()", "mat204", FatalException, ed);
308 }
309
310 if (photonEnergies.size() == 1) {
312 ed << "AddProperty warning. A material property vector must have more than one value.\n"
313 << "Unless you will later add an entry, this is an error.\n"
314 << "Property name: " << key;
315 G4Exception("G4MaterialPropertiesTable::AddProperty()", "mat218", JustWarning, ed);
316 }
317
318 // G4PhysicsVector assumes energies are in increasing order
319 for (std::size_t i = 0; i < photonEnergies.size() - 1; ++i) {
320 if (photonEnergies.at(i + 1) < photonEnergies.at(i)) {
322 ed << "Energies in material property vector must be in increasing "
323 << "order. Key: " << key << " Energy: " << photonEnergies.at(i + 1);
324 G4Exception("G4MaterialPropertiesTable::AddProperty()", "mat215", FatalException, ed);
325 }
326 }
327
328 // if the key doesn't exist, add it if requested
329 if (std::find(fMatPropNames.cbegin(), fMatPropNames.cend(), key) == fMatPropNames.cend()) {
330 if (createNewKey) {
331 fMatPropNames.push_back(key);
332 fMP.push_back(nullptr);
333 }
334 else {
336 ed << "Attempting to create a new material property vector " << key << " without setting\n"
337 << "createNewKey parameter of AddProperty to true.";
338 G4Exception("G4MaterialPropertiesTable::AddProperty()", "mat205", FatalException, ed);
339 }
340 }
341
342 auto* mpv = new G4MaterialPropertyVector(photonEnergies, propertyValues, spline);
343 mpv->SetVerboseLevel(1);
344 if (spline) {
345 mpv->FillSecondDerivatives();
346 }
347 G4int index = GetPropertyIndex(key);
348 fMP[index] = mpv;
349
350 // if key is RINDEX, we calculate GROUPVEL -
351 // contribution from Tao Lin (IHEP, the JUNO experiment)
352 if (key == "RINDEX") {
353 CalculateGROUPVEL();
354 }
355
356 return mpv;
357}
@ JustWarning
G4PhysicsFreeVector G4MaterialPropertyVector

Referenced by AddProperty(), AddProperty(), AddProperty(), G4GDMLReadMaterials::PropertyRead(), G4GDMLReadSolids::PropertyRead(), and G4XrayReflection::SaveHenkeDataAsMaterialProperty().

◆ AddProperty() [5/5]

void G4MaterialPropertiesTable::AddProperty ( const G4String & key,
G4MaterialPropertyVector * opv,
G4bool createNewKey = false )

Definition at line 372 of file G4MaterialPropertiesTable.cc.

374{
375 // Provides a way of adding a property to the Material Properties
376 // Table given an G4MaterialPropertyVector Reference and a key
377
378 // G4PhysicsVector assumes energies are in increasing order
379 if (mpv->GetVectorLength() > 1) {
380 for (std::size_t i = 0; i < mpv->GetVectorLength() - 1; ++i) {
381 if (mpv->Energy(i + 1) < mpv->Energy(i)) {
383 ed << "Energies in material property vector must be in increasing "
384 << "order. Key: " << key << " Energy: " << mpv->Energy(i + 1);
385 G4Exception("G4MaterialPropertiesTable::AddProperty()", "mat216", FatalException, ed);
386 }
387 }
388 }
389
390 if (mpv->GetVectorLength() <= 1) {
392 ed << "AddProperty warning. A material property vector must have more than one value.\n"
393 << "Unless you will later add an entry, this is an error.\n"
394 << "Property name: " << key;
395 G4Exception("G4MaterialPropertiesTable::AddProperty()", "mat219", JustWarning, ed);
396 }
397
398 // if the key doesn't exist, add it
399 if (std::find(fMatPropNames.cbegin(), fMatPropNames.cend(), key) == fMatPropNames.cend()) {
400 if (createNewKey) {
401 fMatPropNames.push_back(key);
402 fMP.push_back(nullptr);
403 }
404 else {
406 ed << "Attempting to create a new material property key " << key << " without setting\n"
407 << "createNewKey parameter of AddProperty to true.";
408 G4Exception("G4MaterialPropertiesTable::AddProperty()", "mat206", FatalException, ed);
409 }
410 }
411 G4int index = GetPropertyIndex(key);
412 fMP[index] = mpv;
413
414 // if key is RINDEX, we calculate GROUPVEL -
415 // contribution from Tao Lin (IHEP, the JUNO experiment)
416 if (key == "RINDEX") {
417 CalculateGROUPVEL();
418 }
419}

◆ ConstPropertyExists() [1/3]

G4bool G4MaterialPropertiesTable::ConstPropertyExists ( const char * key) const

Definition at line 260 of file G4MaterialPropertiesTable.cc.

261{
262 std::size_t index = std::distance(fMatConstPropNames.cbegin(),
263 std::find(fMatConstPropNames.cbegin(), fMatConstPropNames.cend(), key));
264 if (index < fMatConstPropNames.size()) { // index is type std::size_t so >= 0
265 return ConstPropertyExists((G4int)index);
266 }
267 return false;
268}
G4bool ConstPropertyExists(const G4String &key) const

◆ ConstPropertyExists() [2/3]

G4bool G4MaterialPropertiesTable::ConstPropertyExists ( const G4int index) const

Definition at line 242 of file G4MaterialPropertiesTable.cc.

243{
244 // Returns true if a const property corresponding to 'index' exists
245
246 return index >= 0 && index < (G4int)fMCP.size() && fMCP[index].second;
247}

◆ ConstPropertyExists() [3/3]

G4bool G4MaterialPropertiesTable::ConstPropertyExists ( const G4String & key) const

Definition at line 249 of file G4MaterialPropertiesTable.cc.

250{
251 // Returns true if a const property 'key' exists
252 std::size_t index = std::distance(fMatConstPropNames.cbegin(),
253 std::find(fMatConstPropNames.cbegin(), fMatConstPropNames.cend(), key));
254 if (index < fMatConstPropNames.size()) { // index is type std::size_t so >= 0
255 return ConstPropertyExists((G4int)index);
256 }
257 return false;
258}

Referenced by ConstPropertyExists(), ConstPropertyExists(), G4Scintillation::GetScintillationYieldByParticleType(), G4UCNMaterialPropertiesTable::InitMicroRoughnessTables(), G4OpBoundaryProcess::PostStepDoIt(), G4OpWLS2::PostStepDoIt(), G4OpWLS::PostStepDoIt(), G4Scintillation::PostStepDoIt(), and G4UCNMaterialPropertiesTable::SetMicroRoughnessParameters().

◆ DumpTable()

void G4MaterialPropertiesTable::DumpTable ( ) const

Definition at line 531 of file G4MaterialPropertiesTable.cc.

532{
533 // material properties
534 G4int j = 0;
535 for (const auto& prop : fMP) {
536 if (prop != nullptr) {
537 G4cout << j << ": " << fMatPropNames[j] << G4endl;
538 prop->DumpValues();
539 }
540 ++j;
541 }
542 // material constant properties
543 j = 0;
544 for (const auto& cprop : fMCP) {
545 if (cprop.second) {
546 G4cout << j << ": " << fMatConstPropNames[j] << " " << cprop.first << G4endl;
547 }
548 ++j;
549 }
550}
#define G4endl
Definition G4ios.hh:67
G4GLOB_DLL std::ostream G4cout

◆ GetConstProperties()

const std::vector< std::pair< G4double, G4bool > > & G4MaterialPropertiesTable::GetConstProperties ( ) const
inline

Definition at line 141 of file G4MaterialPropertiesTable.hh.

141{ return fMCP; }

Referenced by G4GDMLWriteMaterials::PropertyWrite(), and G4GDMLWriteSolids::PropertyWrite().

◆ GetConstProperty() [1/3]

G4double G4MaterialPropertiesTable::GetConstProperty ( const char * key) const

Definition at line 237 of file G4MaterialPropertiesTable.cc.

238{
240}
G4double GetConstProperty(const G4String &key) const

◆ GetConstProperty() [2/3]

G4double G4MaterialPropertiesTable::GetConstProperty ( const G4int index) const

Definition at line 215 of file G4MaterialPropertiesTable.cc.

216{
217 // Returns the constant material property corresponding to an index
218 // fatal exception if property not found
219
220 if (index < (G4int)fMCP.size() && fMCP[index].second) {
221 return fMCP[index].first;
222 }
224 ed << "Constant Material Property " << fMatConstPropNames[index] << " not found.";
225 G4Exception("G4MaterialPropertiesTable::GetConstProperty()", "mat202", FatalException, ed);
226 return 0.;
227}

◆ GetConstProperty() [3/3]

◆ GetConstPropertyIndex()

G4int G4MaterialPropertiesTable::GetConstPropertyIndex ( const G4String & key) const

Definition at line 185 of file G4MaterialPropertiesTable.cc.

186{
187 // Returns the constant material property index corresponding to a key
188
189 std::size_t index = std::distance(fMatConstPropNames.cbegin(),
190 std::find(fMatConstPropNames.cbegin(), fMatConstPropNames.cend(), key));
191 if (index < fMatConstPropNames.size()) {
192 return (G4int)index;
193 }
194
196 ed << "Constant Material Property Index for key " << key << " not found.";
197 G4Exception("G4MaterialPropertiesTable::GetConstPropertyIndex()", "mat200", FatalException, ed);
198 return 0;
199}

Referenced by AddConstProperty(), GetConstProperty(), GetConstProperty(), and RemoveConstProperty().

◆ GetMaterialConstPropertyNames()

const std::vector< G4String > & G4MaterialPropertiesTable::GetMaterialConstPropertyNames ( ) const
inline

Definition at line 138 of file G4MaterialPropertiesTable.hh.

138{ return fMatConstPropNames; }

Referenced by G4GDMLWriteMaterials::PropertyWrite(), and G4GDMLWriteSolids::PropertyWrite().

◆ GetMaterialPropertyNames()

const std::vector< G4String > & G4MaterialPropertiesTable::GetMaterialPropertyNames ( ) const
inline

Definition at line 137 of file G4MaterialPropertiesTable.hh.

137{ return fMatPropNames; }

Referenced by G4GDMLWriteMaterials::PropertyWrite(), and G4GDMLWriteSolids::PropertyWrite().

◆ GetProperties()

const std::vector< G4MaterialPropertyVector * > & G4MaterialPropertiesTable::GetProperties ( ) const
inline

Definition at line 140 of file G4MaterialPropertiesTable.hh.

140{ return fMP; }

Referenced by G4GDMLWriteMaterials::PropertyWrite(), and G4GDMLWriteSolids::PropertyWrite().

◆ GetProperty() [1/3]

◆ GetProperty() [2/3]

G4MaterialPropertyVector * G4MaterialPropertiesTable::GetProperty ( const G4int index) const

Definition at line 289 of file G4MaterialPropertiesTable.cc.

290{
291 // Returns a Material Property Vector corresponding to an index
292 // returns nullptr if the property has not been defined by user
293 if (index >= 0 && index < (G4int)fMP.size()) {
294 return fMP[index];
295 }
296 return nullptr;
297}

◆ GetProperty() [3/3]

G4MaterialPropertyVector * G4MaterialPropertiesTable::GetProperty ( const G4String & key) const

Definition at line 270 of file G4MaterialPropertiesTable.cc.

271{
272 // Returns a Material Property Vector corresponding to a key
273 if (std::find(fMatPropNames.cbegin(), fMatPropNames.cend(), key) != fMatPropNames.cend()) {
274 const G4int index = GetPropertyIndex(G4String(key));
275 return GetProperty(index);
276 }
277 return nullptr;
278}

◆ GetPropertyIndex()

G4int G4MaterialPropertiesTable::GetPropertyIndex ( const G4String & key) const

Definition at line 201 of file G4MaterialPropertiesTable.cc.

202{
203 // Returns the material property index corresponding to a key
204 std::size_t index = std::distance(
205 fMatPropNames.cbegin(), std::find(fMatPropNames.cbegin(), fMatPropNames.cend(), key));
206 if (index < fMatPropNames.size()) {
207 return (G4int)index;
208 }
210 ed << "Material Property Index for key " << key << " not found.";
211 G4Exception("G4MaterialPropertiesTable::GetPropertyIndex()", "mat201", FatalException, ed);
212 return 0;
213}

Referenced by AddEntry(), AddProperty(), AddProperty(), GetProperty(), GetProperty(), and RemoveProperty().

◆ RemoveConstProperty() [1/2]

void G4MaterialPropertiesTable::RemoveConstProperty ( const char * key)

Definition at line 475 of file G4MaterialPropertiesTable.cc.

476{
478}
void RemoveConstProperty(const G4String &key)

◆ RemoveConstProperty() [2/2]

void G4MaterialPropertiesTable::RemoveConstProperty ( const G4String & key)

Definition at line 467 of file G4MaterialPropertiesTable.cc.

468{
469 G4int index = GetConstPropertyIndex(key);
470 if (index < (G4int)fMCP.size()) {
471 fMCP[index] = std::pair<G4double, G4bool>{0., false};
472 }
473}

Referenced by RemoveConstProperty(), and G4UCNMaterialPropertiesTable::SetMicroRoughnessParameters().

◆ RemoveProperty() [1/2]

void G4MaterialPropertiesTable::RemoveProperty ( const char * key)

Definition at line 487 of file G4MaterialPropertiesTable.cc.

487{ RemoveProperty(G4String(key)); }
void RemoveProperty(const G4String &key)

◆ RemoveProperty() [2/2]

void G4MaterialPropertiesTable::RemoveProperty ( const G4String & key)

Definition at line 480 of file G4MaterialPropertiesTable.cc.

481{
482 G4int index = GetPropertyIndex(key);
483 delete fMP[index];
484 fMP[index] = nullptr;
485}

Referenced by RemoveProperty().


The documentation for this class was generated from the following files: