| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 << 25 //
>> 26 //
>> 27 // $Id: G4MaterialPropertiesTable.cc,v 1.19 2007/04/25 15:46:55 gum Exp $
>> 28 // GEANT4 tag $Name: geant4-09-01-patch-01 $
>> 29 //
>> 30 //
26 ////////////////////////////////////////////// 31 ////////////////////////////////////////////////////////////////////////
27 // G4MaterialPropertiesTable Implementation 32 // G4MaterialPropertiesTable Implementation
28 ////////////////////////////////////////////// 33 ////////////////////////////////////////////////////////////////////////
29 // 34 //
30 // File: G4MaterialPropertiesTable.cc << 35 // File: G4MaterialPropertiesTable.cc
31 // Version: 1.0 36 // Version: 1.0
32 // Created: 1996-02-08 37 // Created: 1996-02-08
33 // Author: Juliet Armstrong 38 // Author: Juliet Armstrong
34 // Updated: 2005-05-12 add SetGROUPVEL(), 39 // Updated: 2005-05-12 add SetGROUPVEL(), courtesy of
35 // Horton-Smith (bug report #741) 40 // Horton-Smith (bug report #741), by P. Gumplinger
36 // 2002-11-05 add named material 41 // 2002-11-05 add named material constants by P. Gumplinger
37 // 1999-11-05 Migration from G4RW 42 // 1999-11-05 Migration from G4RWTPtrHashDictionary to STL
38 // by John Allison 43 // by John Allison
39 // 1997-03-26 by Peter Gumplinger 44 // 1997-03-26 by Peter Gumplinger
40 // > cosmetics (only) 45 // > cosmetics (only)
>> 46 // mail: gum@triumf.ca
41 // 47 //
42 ////////////////////////////////////////////// 48 ////////////////////////////////////////////////////////////////////////
43 49
44 #include "G4MaterialPropertiesTable.hh" <<
45 <<
46 #include "G4Log.hh" <<
47 #include "G4OpticalMaterialProperties.hh" <<
48 #include "G4PhysicalConstants.hh" <<
49 #include "globals.hh" 50 #include "globals.hh"
>> 51 #include "G4MaterialPropertiesTable.hh"
50 52
51 #include <algorithm> << 53 /////////////////
52 #include <cassert> << 54 // Constructors
53 << 55 /////////////////
54 #ifdef G4MULTITHREADED << 56
55 # include "G4AutoLock.hh" << 57 G4MaterialPropertiesTable::G4MaterialPropertiesTable() {}
56 namespace << 58
57 { << 59 ////////////////
58 G4Mutex materialPropertyTableMutex = G4MUTEX_I << 60 // Destructors
59 } << 61 ////////////////
60 #endif <<
61 <<
62 G4MaterialPropertiesTable::G4MaterialPropertie <<
63 { <<
64 // elements of these 2 vectors must be in sa <<
65 // the corresponding enums in G4MaterialProp <<
66 fMatPropNames.assign(kNumberOfPropertyIndex, <<
67 fMatPropNames[kRINDEX] = <<
68 fMatPropNames[kREFLECTIVITY] = <<
69 fMatPropNames[kREALRINDEX] = <<
70 fMatPropNames[kIMAGINARYRINDEX] = <<
71 fMatPropNames[kEFFICIENCY] = <<
72 fMatPropNames[kTRANSMITTANCE] = <<
73 fMatPropNames[kSPECULARLOBECONSTANT] = <<
74 fMatPropNames[kSPECULARSPIKECONSTANT] = <<
75 fMatPropNames[kBACKSCATTERCONSTANT] = <<
76 fMatPropNames[kGROUPVEL] = <<
77 fMatPropNames[kMIEHG] = <<
78 fMatPropNames[kRAYLEIGH] = <<
79 fMatPropNames[kWLSCOMPONENT] = <<
80 fMatPropNames[kWLSABSLENGTH] = <<
81 fMatPropNames[kWLSCOMPONENT2] = <<
82 fMatPropNames[kWLSABSLENGTH2] = <<
83 fMatPropNames[kABSLENGTH] = <<
84 fMatPropNames[kPROTONSCINTILLATIONYIELD] = <<
85 fMatPropNames[kDEUTERONSCINTILLATIONYIELD] = <<
86 fMatPropNames[kTRITONSCINTILLATIONYIELD] = <<
87 fMatPropNames[kALPHASCINTILLATIONYIELD] = <<
88 fMatPropNames[kIONSCINTILLATIONYIELD] = <<
89 fMatPropNames[kELECTRONSCINTILLATIONYIELD] = <<
90 fMatPropNames[kSCINTILLATIONCOMPONENT1] = <<
91 fMatPropNames[kSCINTILLATIONCOMPONENT2] = <<
92 fMatPropNames[kSCINTILLATIONCOMPONENT3] = <<
93 fMatPropNames[kCOATEDRINDEX] = <<
94 <<
95 fMP.assign(kNumberOfPropertyIndex, nullptr); <<
96 <<
97 fMatConstPropNames.assign(kNumberOfConstProp <<
98 fMatConstPropNames[kSURFACEROUGHNESS] = <<
99 fMatConstPropNames[kISOTHERMAL_COMPRESSIBILI <<
100 fMatConstPropNames[kRS_SCALE_FACTOR] = <<
101 fMatConstPropNames[kWLSMEANNUMBERPHOTONS] = <<
102 fMatConstPropNames[kWLSTIMECONSTANT] = <<
103 fMatConstPropNames[kWLSMEANNUMBERPHOTONS2] = <<
104 fMatConstPropNames[kWLSTIMECONSTANT2] = <<
105 fMatConstPropNames[kMIEHG_FORWARD] = <<
106 fMatConstPropNames[kMIEHG_BACKWARD] = <<
107 fMatConstPropNames[kMIEHG_FORWARD_RATIO] = <<
108 fMatConstPropNames[kSCINTILLATIONYIELD] = <<
109 fMatConstPropNames[kRESOLUTIONSCALE] = <<
110 fMatConstPropNames[kFERMIPOT] = <<
111 fMatConstPropNames[kDIFFUSION] = <<
112 fMatConstPropNames[kSPINFLIP] = <<
113 fMatConstPropNames[kLOSS] = <<
114 fMatConstPropNames[kLOSSCS] = <<
115 fMatConstPropNames[kABSCS] = <<
116 fMatConstPropNames[kSCATCS] = <<
117 fMatConstPropNames[kMR_NBTHETA] = <<
118 fMatConstPropNames[kMR_NBE] = <<
119 fMatConstPropNames[kMR_RRMS] = <<
120 fMatConstPropNames[kMR_CORRLEN] = <<
121 fMatConstPropNames[kMR_THETAMIN] = <<
122 fMatConstPropNames[kMR_THETAMAX] = <<
123 fMatConstPropNames[kMR_EMIN] = <<
124 fMatConstPropNames[kMR_EMAX] = <<
125 fMatConstPropNames[kMR_ANGNOTHETA] = <<
126 fMatConstPropNames[kMR_ANGNOPHI] = <<
127 fMatConstPropNames[kMR_ANGCUT] = <<
128 fMatConstPropNames[kSCINTILLATIONTIMECONSTAN <<
129 fMatConstPropNames[kSCINTILLATIONTIMECONSTAN <<
130 fMatConstPropNames[kSCINTILLATIONTIMECONSTAN <<
131 fMatConstPropNames[kSCINTILLATIONRISETIME1] <<
132 fMatConstPropNames[kSCINTILLATIONRISETIME2] <<
133 fMatConstPropNames[kSCINTILLATIONRISETIME3] <<
134 fMatConstPropNames[kSCINTILLATIONYIELD1] = <<
135 fMatConstPropNames[kSCINTILLATIONYIELD2] = <<
136 fMatConstPropNames[kSCINTILLATIONYIELD3] = <<
137 fMatConstPropNames[kPROTONSCINTILLATIONYIELD <<
138 fMatConstPropNames[kPROTONSCINTILLATIONYIELD <<
139 fMatConstPropNames[kPROTONSCINTILLATIONYIELD <<
140 fMatConstPropNames[kDEUTERONSCINTILLATIONYIE <<
141 fMatConstPropNames[kDEUTERONSCINTILLATIONYIE <<
142 fMatConstPropNames[kDEUTERONSCINTILLATIONYIE <<
143 fMatConstPropNames[kTRITONSCINTILLATIONYIELD <<
144 fMatConstPropNames[kTRITONSCINTILLATIONYIELD <<
145 fMatConstPropNames[kTRITONSCINTILLATIONYIELD <<
146 fMatConstPropNames[kALPHASCINTILLATIONYIELD1 <<
147 fMatConstPropNames[kALPHASCINTILLATIONYIELD2 <<
148 fMatConstPropNames[kALPHASCINTILLATIONYIELD3 <<
149 fMatConstPropNames[kIONSCINTILLATIONYIELD1] <<
150 fMatConstPropNames[kIONSCINTILLATIONYIELD2] <<
151 fMatConstPropNames[kIONSCINTILLATIONYIELD3] <<
152 fMatConstPropNames[kELECTRONSCINTILLATIONYIE <<
153 fMatConstPropNames[kELECTRONSCINTILLATIONYIE <<
154 fMatConstPropNames[kELECTRONSCINTILLATIONYIE <<
155 fMatConstPropNames[kCOATEDTHICKNESS] = <<
156 fMatConstPropNames[kCOATEDFRUSTRATEDTRANSMIS <<
157 fMatConstPropNames[kPROTONSCINTILLATIONTIMEC <<
158 fMatConstPropNames[kPROTONSCINTILLATIONTIMEC <<
159 fMatConstPropNames[kPROTONSCINTILLATIONTIMEC <<
160 fMatConstPropNames[kDEUTERONSCINTILLATIONTIM <<
161 fMatConstPropNames[kDEUTERONSCINTILLATIONTIM <<
162 fMatConstPropNames[kDEUTERONSCINTILLATIONTIM <<
163 fMatConstPropNames[kTRITONSCINTILLATIONTIMEC <<
164 fMatConstPropNames[kTRITONSCINTILLATIONTIMEC <<
165 fMatConstPropNames[kTRITONSCINTILLATIONTIMEC <<
166 fMatConstPropNames[kALPHASCINTILLATIONTIMECO <<
167 fMatConstPropNames[kALPHASCINTILLATIONTIMECO <<
168 fMatConstPropNames[kALPHASCINTILLATIONTIMECO <<
169 fMatConstPropNames[kIONSCINTILLATIONTIMECONS <<
170 fMatConstPropNames[kIONSCINTILLATIONTIMECONS <<
171 fMatConstPropNames[kIONSCINTILLATIONTIMECONS <<
172 fMatConstPropNames[kELECTRONSCINTILLATIONTIM <<
173 fMatConstPropNames[kELECTRONSCINTILLATIONTIM <<
174 fMatConstPropNames[kELECTRONSCINTILLATIONTIM <<
175 fMCP.assign(kNumberOfConstPropertyIndex, {0. <<
176 } <<
177 62
178 G4MaterialPropertiesTable::~G4MaterialProperti 63 G4MaterialPropertiesTable::~G4MaterialPropertiesTable()
179 { 64 {
180 for (auto prop : fMP) { << 65 MPTiterator i;
181 delete (prop); << 66 for (i = MPT.begin(); i != MPT.end(); ++i) {
182 } << 67 delete (*i).second;
>> 68 }
>> 69 MPT.clear();
>> 70 MPTC.clear();
183 } 71 }
184 72
185 G4int G4MaterialPropertiesTable::GetConstPrope << 73 ////////////
186 { << 74 // Methods
187 // Returns the constant material property in << 75 ////////////
188 <<
189 std::size_t index = std::distance(fMatConstP <<
190 std::find(fMatConstPropNames.cbegin(), fMa <<
191 if (index < fMatConstPropNames.size()) { <<
192 return (G4int)index; <<
193 } <<
194 <<
195 G4ExceptionDescription ed; <<
196 ed << "Constant Material Property Index for <<
197 G4Exception("G4MaterialPropertiesTable::GetC <<
198 return 0; <<
199 } <<
200 76
201 G4int G4MaterialPropertiesTable::GetPropertyIn << 77 void G4MaterialPropertiesTable::AddConstProperty(const char *key,
202 { << 78 G4double PropertyValue)
203 // Returns the material property index corre <<
204 std::size_t index = std::distance( <<
205 fMatPropNames.cbegin(), std::find(fMatProp <<
206 if (index < fMatPropNames.size()) { <<
207 return (G4int)index; <<
208 } <<
209 G4ExceptionDescription ed; <<
210 ed << "Material Property Index for key " << <<
211 G4Exception("G4MaterialPropertiesTable::GetP <<
212 return 0; <<
213 } <<
214 <<
215 G4double G4MaterialPropertiesTable::GetConstPr <<
216 { <<
217 // Returns the constant material property co <<
218 // fatal exception if property not found <<
219 <<
220 if (index < (G4int)fMCP.size() && fMCP[index <<
221 return fMCP[index].first; <<
222 } <<
223 G4ExceptionDescription ed; <<
224 ed << "Constant Material Property " << fMatC <<
225 G4Exception("G4MaterialPropertiesTable::GetC <<
226 return 0.; <<
227 } <<
228 <<
229 G4double G4MaterialPropertiesTable::GetConstPr <<
230 { 79 {
231 // Returns the constant material property co << 80 // Provides a way of adding a constant property to the Mataerial Properties
232 // fatal exception if property not found << 81 // Table given a key
233 82
234 return GetConstProperty(GetConstPropertyInde << 83 MPTC [G4String(key)] = PropertyValue;
235 } <<
236 84
237 G4double G4MaterialPropertiesTable::GetConstPr <<
238 { <<
239 return GetConstProperty(GetConstPropertyInde <<
240 } 85 }
241 86
242 G4bool G4MaterialPropertiesTable::ConstPropert << 87 void G4MaterialPropertiesTable::AddProperty(const char *key,
>> 88 G4double *PhotonMomenta,
>> 89 G4double *PropertyValues,
>> 90 G4int NumEntries)
243 { 91 {
244 // Returns true if a const property correspo << 92 // Privides a way of adding a property to the Material Properties
>> 93 // Table given a pair of numbers and a key
245 94
246 return index >= 0 && index < (G4int)fMCP.siz << 95 G4MaterialPropertyVector *mpv =
247 } << 96 new G4MaterialPropertyVector(PhotonMomenta,
248 << 97 PropertyValues,
249 G4bool G4MaterialPropertiesTable::ConstPropert << 98 NumEntries);
250 { << 99 MPT [G4String(key)] = mpv;
251 // Returns true if a const property 'key' ex <<
252 std::size_t index = std::distance(fMatConstP <<
253 std::find(fMatConstPropNames.cbegin(), fMa <<
254 if (index < fMatConstPropNames.size()) { // <<
255 return ConstPropertyExists((G4int)index); <<
256 } <<
257 return false; <<
258 } 100 }
259 101
260 G4bool G4MaterialPropertiesTable::ConstPropert << 102 void G4MaterialPropertiesTable::AddProperty(const char *key,
>> 103 G4MaterialPropertyVector *mpv)
261 { 104 {
262 std::size_t index = std::distance(fMatConstP << 105 // Provides a way of adding a property to the Material Properties
263 std::find(fMatConstPropNames.cbegin(), fMa << 106 // Table given an G4MaterialPropertyVector Reference and a key
264 if (index < fMatConstPropNames.size()) { // <<
265 return ConstPropertyExists((G4int)index); <<
266 } <<
267 return false; <<
268 } <<
269 107
270 G4MaterialPropertyVector* G4MaterialProperties << 108 MPT [G4String(key)] = mpv;
271 { << 109 }
272 // Returns a Material Property Vector corres <<
273 if (std::find(fMatPropNames.cbegin(), fMatPr <<
274 const G4int index = GetPropertyIndex(G4Str <<
275 return GetProperty(index); <<
276 } <<
277 return nullptr; <<
278 } <<
279 110
280 G4MaterialPropertyVector* G4MaterialProperties << 111 void G4MaterialPropertiesTable::RemoveConstProperty(const char *key)
281 { 112 {
282 if (std::find(fMatPropNames.cbegin(), fMatPr << 113 MPTC.erase(G4String(key));
283 const G4int index = GetPropertyIndex(G4Str <<
284 return GetProperty(index); <<
285 } <<
286 return nullptr; <<
287 } 114 }
288 115
289 G4MaterialPropertyVector* G4MaterialProperties << 116 void G4MaterialPropertiesTable::RemoveProperty(const char *key)
290 { 117 {
291 // Returns a Material Property Vector corres << 118 MPT.erase(G4String(key));
292 // returns nullptr if the property has not b <<
293 if (index >= 0 && index < (G4int)fMP.size()) <<
294 return fMP[index]; <<
295 } <<
296 return nullptr; <<
297 } 119 }
298 120
299 G4MaterialPropertyVector* G4MaterialProperties << 121 G4double G4MaterialPropertiesTable::GetConstProperty(const char *key)
300 const std::vector<G4double>& photonEnergies, <<
301 G4bool createNewKey, G4bool spline) <<
302 { <<
303 if (photonEnergies.size() != propertyValues. <<
304 G4ExceptionDescription ed; <<
305 ed << "AddProperty error. Number of proper <<
306 << "energy values. Property name: " << <<
307 G4Exception("G4MaterialPropertiesTable::Ad <<
308 } <<
309 <<
310 if (photonEnergies.size() == 1) { <<
311 G4ExceptionDescription ed; <<
312 ed << "AddProperty warning. A material pro <<
313 << "Unless you will later add an entry, <<
314 << "Property name: " << key; <<
315 G4Exception("G4MaterialPropertiesTable::Ad <<
316 } <<
317 <<
318 // G4PhysicsVector assumes energies are in i <<
319 for (std::size_t i = 0; i < photonEnergies.s <<
320 if (photonEnergies.at(i + 1) < photonEnerg <<
321 G4ExceptionDescription ed; <<
322 ed << "Energies in material property vec <<
323 << "order. Key: " << key << " Energy: <<
324 G4Exception("G4MaterialPropertiesTable:: <<
325 } <<
326 } <<
327 <<
328 // if the key doesn't exist, add it if reque <<
329 if (std::find(fMatPropNames.cbegin(), fMatPr <<
330 if (createNewKey) { <<
331 fMatPropNames.push_back(key); <<
332 fMP.push_back(nullptr); <<
333 } <<
334 else { <<
335 G4ExceptionDescription ed; <<
336 ed << "Attempting to create a new materi <<
337 << "createNewKey parameter of AddProp <<
338 G4Exception("G4MaterialPropertiesTable:: <<
339 } <<
340 } <<
341 <<
342 auto* mpv = new G4MaterialPropertyVector(pho <<
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 <<
352 if (key == "RINDEX") { <<
353 CalculateGROUPVEL(); <<
354 } <<
355 <<
356 return mpv; <<
357 } <<
358 <<
359 G4MaterialPropertyVector* G4MaterialProperties <<
360 G4double* photonEnergies, G4double* property <<
361 G4bool spline) <<
362 { 122 {
363 // Provides a way of adding a property to th << 123 // Returns the constant material property corresponding to a key
364 // Table given a pair of arrays and a key <<
365 G4String k(key); <<
366 124
367 std::vector<G4double> energies(photonEnergie << 125 MPTCiterator j;
368 std::vector<G4double> values(propertyValues, << 126 j = MPTC.find(G4String(key));
369 return AddProperty(k, energies, values, crea << 127 if ( j != MPTC.end() ) {
>> 128 return j->second;
>> 129 }
>> 130 else {
>> 131 G4Exception("G4MaterialPropertiesTable::GetConstProperty ==> "
>> 132 "Constant Material Property not found.");
>> 133 return G4double(0.0);
>> 134 }
370 } 135 }
371 136
372 void G4MaterialPropertiesTable::AddProperty( << 137 G4bool G4MaterialPropertiesTable::ConstPropertyExists(const char *key)
373 const G4String& key, G4MaterialPropertyVecto <<
374 { 138 {
375 // Provides a way of adding a property to t << 139 // Return true if a const property 'key' exists
376 // Table given an G4MaterialPropertyVector <<
377 <<
378 // G4PhysicsVector assumes energies are in i <<
379 if (mpv->GetVectorLength() > 1) { <<
380 for (std::size_t i = 0; i < mpv->GetVector <<
381 if (mpv->Energy(i + 1) < mpv->Energy(i)) <<
382 G4ExceptionDescription ed; <<
383 ed << "Energies in material property v <<
384 << "order. Key: " << key << " Energ <<
385 G4Exception("G4MaterialPropertiesTable <<
386 } <<
387 } <<
388 } <<
389 <<
390 if (mpv->GetVectorLength() <= 1) { <<
391 G4ExceptionDescription ed; <<
392 ed << "AddProperty warning. A material pro <<
393 << "Unless you will later add an entry, <<
394 << "Property name: " << key; <<
395 G4Exception("G4MaterialPropertiesTable::Ad <<
396 } <<
397 140
398 // if the key doesn't exist, add it << 141 MPTCiterator j;
399 if (std::find(fMatPropNames.cbegin(), fMatPr << 142 j = MPTC.find(G4String(key));
400 if (createNewKey) { << 143 if ( j != MPTC.end() ) {
401 fMatPropNames.push_back(key); << 144 return true;
402 fMP.push_back(nullptr); << 145 }
403 } << 146 else {
404 else { << 147 return false;
405 G4ExceptionDescription ed; << 148 }
406 ed << "Attempting to create a new materi <<
407 << "createNewKey parameter of AddProp <<
408 G4Exception("G4MaterialPropertiesTable:: <<
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 <<
416 if (key == "RINDEX") { <<
417 CalculateGROUPVEL(); <<
418 } <<
419 } 149 }
420 150
421 void G4MaterialPropertiesTable::AddProperty( << 151 G4MaterialPropertyVector* G4MaterialPropertiesTable::GetProperty(const char *key)
422 const char* key, G4MaterialPropertyVector* m <<
423 { 152 {
424 AddProperty(G4String(key), mpv, createNewKey << 153 // Returns a Material Property Vector corresponding to a key
425 } <<
426 154
427 void G4MaterialPropertiesTable::AddProperty(co << 155 if(MPT[G4String(key)] != 0) {
428 { << 156 return MPT[G4String(key)];
429 // load a material property vector defined i << 157 }else{
430 G4MaterialPropertyVector* v = G4OpticalMater << 158 if(G4String(key) == "GROUPVEL") {
431 AddProperty(key, v); << 159 return SetGROUPVEL();
>> 160 }else{
>> 161 return MPT[G4String(key)];
>> 162 }
>> 163 }
432 } 164 }
433 165
434 void G4MaterialPropertiesTable::AddConstProper << 166 void G4MaterialPropertiesTable::AddEntry(const char *key,
435 const G4String& key, G4double propertyValue, << 167 G4double aPhotonMomentum,
436 { << 168 G4double aPropertyValue)
437 // Provides a way of adding a constant prope <<
438 // Table given a key <<
439 if (std::find(fMatConstPropNames.cbegin(), f <<
440 fMatConstPropNames.cend()) <<
441 { <<
442 if (createNewKey) { <<
443 fMatConstPropNames.push_back(key); <<
444 fMCP.emplace_back(0., true); <<
445 } <<
446 else { <<
447 G4ExceptionDescription ed; <<
448 ed << "Attempting to create a new materi <<
449 << " without setting\n" <<
450 << "createNewKey parameter of AddProp <<
451 G4Exception("G4MaterialPropertiesTable:: <<
452 } <<
453 } <<
454 G4int index = GetConstPropertyIndex(key); <<
455 169
456 fMCP[index] = std::pair<G4double, G4bool>{pr << 170 // Allows to add an entry pair directly to the Material Property Vector
457 } << 171 // given a key
458 172
459 void G4MaterialPropertiesTable::AddConstProper <<
460 const char* key, G4double propertyValue, G4b <<
461 { 173 {
462 // Provides a way of adding a constant prope << 174 G4MaterialPropertyVector *targetVector=MPT [G4String(key)];
463 // Table given a key << 175 if (targetVector != 0) {
464 AddConstProperty(G4String(key), propertyValu << 176 targetVector->AddElement(aPhotonMomentum, aPropertyValue);
>> 177 }
>> 178 else {
>> 179 G4Exception("G4MaterialPropertiesTable::AddEntry ==> "
>> 180 "Material Property Vector not found.");
>> 181 }
465 } 182 }
466 183
467 void G4MaterialPropertiesTable::RemoveConstPro << 184 void G4MaterialPropertiesTable::RemoveEntry(const char *key,
>> 185 G4double aPhotonMomentum)
468 { 186 {
469 G4int index = GetConstPropertyIndex(key); << 187 // Allows to remove an entry pair directly from the Material Property Vector
470 if (index < (G4int)fMCP.size()) { << 188 // given a key
471 fMCP[index] = std::pair<G4double, G4bool>{ <<
472 } <<
473 } <<
474 189
475 void G4MaterialPropertiesTable::RemoveConstPro << 190 G4MaterialPropertyVector *targetVector=MPT [G4String(key)];
476 { << 191 if (targetVector) {
477 RemoveConstProperty(G4String(key)); << 192 targetVector->RemoveElement(aPhotonMomentum);
>> 193 }
>> 194 else {
>> 195 G4Exception("G4MaterialPropertiesTable::RemoveEntry ==> "
>> 196 "Material Property Vector not found.");
>> 197 }
478 } 198 }
479 << 199 void G4MaterialPropertiesTable::DumpTable()
480 void G4MaterialPropertiesTable::RemoveProperty <<
481 { 200 {
482 G4int index = GetPropertyIndex(key); << 201 MPTiterator i;
483 delete fMP[index]; << 202 for (i = MPT.begin(); i != MPT.end(); ++i) {
484 fMP[index] = nullptr; << 203 G4cout << (*i).first << G4endl;
485 } << 204 if ( (*i).second != 0 ) {
486 << 205 (*i).second->DumpVector();
487 void G4MaterialPropertiesTable::RemoveProperty << 206 }
488 << 207 else {
489 void G4MaterialPropertiesTable::AddEntry( << 208 G4cout << "NULL Material Property Vector Pointer." << G4endl;
490 const G4String& key, G4double aPhotonEnergy, << 209 }
491 { <<
492 // Allows to add an entry pair directly to t <<
493 // given a key. <<
494 if (std::find(fMatPropNames.cbegin(), fMatPr <<
495 G4ExceptionDescription ed; <<
496 ed << "Material Property Vector " << key < <<
497 G4Exception("G4MaterialPropertiesTable::Ad <<
498 } <<
499 G4int index = GetPropertyIndex(key); <<
500 <<
501 G4MaterialPropertyVector* targetVector = fMP <<
502 if (targetVector != nullptr) { <<
503 // do not allow duplicate energies <<
504 for (std::size_t i = 0; i < targetVector-> <<
505 if (aPhotonEnergy == targetVector->Energ <<
506 G4ExceptionDescription ed; <<
507 ed << "Energy values in material prope <<
508 << "Key: " << key; <<
509 G4Exception("G4MaterialPropertiesTable <<
510 } <<
511 } <<
512 <<
513 targetVector->InsertValues(aPhotonEnergy, <<
514 } 210 }
515 else { << 211 MPTCiterator j;
516 G4ExceptionDescription ed; << 212 for (j = MPTC.begin(); j != MPTC.end(); ++j) {
517 ed << "Material Property Vector " << key < << 213 G4cout << j->first << G4endl;
518 G4Exception("G4MaterialPropertiesTable::Ad << 214 if ( j->second != 0 ) {
>> 215 G4cout << j->second << G4endl;
>> 216 }
>> 217 else {
>> 218 G4cout << "No Material Constant Property." << G4endl;
>> 219 }
519 } 220 }
520 if (key == "RINDEX") { <<
521 CalculateGROUPVEL(); <<
522 } <<
523 } <<
524 <<
525 void G4MaterialPropertiesTable::AddEntry( <<
526 const char* key, G4double aPhotonEnergy, G4d <<
527 { <<
528 AddEntry(G4String(key), aPhotonEnergy, aProp <<
529 } <<
530 221
531 void G4MaterialPropertiesTable::DumpTable() co <<
532 { <<
533 // material properties <<
534 G4int j = 0; <<
535 for (const auto& prop : fMP) { <<
536 if (prop != nullptr) { <<
537 G4cout << j << ": " << fMatPropNames[j] <<
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 << ": " << fMatConstPropName <<
547 } <<
548 ++j; <<
549 } <<
550 } 222 }
551 << 223 G4MaterialPropertyVector* G4MaterialPropertiesTable::SetGROUPVEL()
552 G4MaterialPropertyVector* G4MaterialProperties <<
553 { 224 {
554 #ifdef G4MULTITHREADED <<
555 G4AutoLock mptm(&materialPropertyTableMutex) <<
556 #endif <<
557 <<
558 // check if "GROUPVEL" already exists. If so <<
559 if (fMP[kGROUPVEL] != nullptr) { <<
560 this->RemoveProperty("GROUPVEL"); <<
561 } <<
562 <<
563 // fetch RINDEX data, give up if unavailable 225 // fetch RINDEX data, give up if unavailable
564 G4MaterialPropertyVector* rindex = this->Get << 226 G4MaterialPropertyVector *rindex = this->GetProperty("RINDEX");
565 if (rindex == nullptr) { << 227 if (rindex==0) return NULL;
566 return nullptr; << 228
567 } << 229 rindex->ResetIterator();
568 <<
569 // RINDEX exists but has no entries, give up 230 // RINDEX exists but has no entries, give up
570 if (rindex->GetVectorLength() == 0) { << 231 if ( (++*rindex) == false ) return NULL;
571 return nullptr; <<
572 } <<
573 232
574 // add GROUPVEL vector 233 // add GROUPVEL vector
575 auto* groupvel = new G4MaterialPropertyVecto << 234 G4MaterialPropertyVector* groupvel = new G4MaterialPropertyVector();
576 groupvel->SetVerboseLevel(1); <<
577 235
>> 236 this->AddProperty( "GROUPVEL", groupvel );
>> 237
578 // fill GROUPVEL vector using RINDEX values 238 // fill GROUPVEL vector using RINDEX values
579 // rindex built-in "iterator" was advanced t 239 // rindex built-in "iterator" was advanced to first entry above
580 G4double E0 = rindex->Energy(0); << 240 G4double E0 = rindex->GetPhotonMomentum();
581 G4double n0 = (*rindex)[0]; << 241 G4double n0 = rindex->GetProperty();
582 242
583 if (E0 <= 0.) { << 243 if (E0 <= 0. ) G4Exception("G4MaterialPropertiesTable::SetGROUPVEL ==> "
584 G4Exception("G4MaterialPropertiesTable::Ca << 244 "Optical Photon Energy <= 0");
585 "Optical Photon Energy <= 0"); << 245
586 } << 246 if ( ++*rindex ) {
587 <<
588 if (rindex->GetVectorLength() >= 2) { <<
589 // good, we have at least two entries in R 247 // good, we have at least two entries in RINDEX
590 // get next energy/value pair 248 // get next energy/value pair
>> 249 G4double E1 = rindex->GetPhotonMomentum();
>> 250 G4double n1 = rindex->GetProperty();
591 251
592 G4double E1 = rindex->Energy(1); << 252 if (E1 <= 0. ) G4Exception("G4MaterialPropertiesTable::SetGROUPVEL ==> "
593 G4double n1 = (*rindex)[1]; << 253 "Optical Photon Energy <= 0");
594 <<
595 if (E1 <= 0.) { <<
596 G4Exception("G4MaterialPropertiesTable:: <<
597 "Optical Photon Energy <= 0"); <<
598 } <<
599 <<
600 G4double vg; 254 G4double vg;
601 <<
602 // add entry at first photon energy 255 // add entry at first photon energy
603 vg = c_light / (n0 + (n1 - n0) / G4Log(E1 << 256 vg = c_light/(n0+(n1-n0)/std::log(E1/E0));
604 <<
605 // allow only for 'normal dispersion' -> d 257 // allow only for 'normal dispersion' -> dn/d(logE) > 0
606 if ((vg < 0) || (vg > c_light / n0)) { << 258 if(vg<0 || vg>c_light/n0)vg = c_light/n0;
607 vg = c_light / n0; << 259 groupvel->AddElement( E0, vg );
608 } <<
609 <<
610 groupvel->InsertValues(E0, vg); <<
611 <<
612 // add entries at midpoints between remain 260 // add entries at midpoints between remaining photon energies
613 for (std::size_t i = 2; i < rindex->GetVec << 261 while(1) {
614 vg = c_light / (0.5 * (n0 + n1) + (n1 - << 262 vg = c_light/( 0.5*(n0+n1)+(n1-n0)/std::log(E1/E0));
615 <<
616 // allow only for 'normal dispersion' -> 263 // allow only for 'normal dispersion' -> dn/d(logE) > 0
617 if ((vg < 0) || (vg > c_light / (0.5 * ( << 264 if(vg<0 || vg>c_light/(0.5*(n0+n1)))vg = c_light/(0.5*(n0+n1));
618 vg = c_light / (0.5 * (n0 + n1)); << 265 groupvel->AddElement( 0.5*(E0+E1), vg );
619 } <<
620 groupvel->InsertValues(0.5 * (E0 + E1), <<
621 <<
622 // get next energy/value pair, or exit l 266 // get next energy/value pair, or exit loop
>> 267 if (!(++*rindex)) break;
623 E0 = E1; 268 E0 = E1;
624 n0 = n1; 269 n0 = n1;
625 E1 = rindex->Energy(i); << 270 E1 = rindex->GetPhotonMomentum();
626 n1 = (*rindex)[i]; << 271 n1 = rindex->GetProperty();
627 272
628 if (E1 <= 0.) { << 273 if (E1 <= 0. ) G4Exception("G4MaterialPropertiesTable::SetGROUPVEL ==> "
629 G4Exception("G4MaterialPropertiesTable << 274 "Optical Photon Energy <= 0");
630 "Optical Photon Energy <= 0"); <<
631 } <<
632 } 275 }
633 <<
634 // add entry at last photon energy 276 // add entry at last photon energy
635 vg = c_light / (n1 + (n1 - n0) / G4Log(E1 << 277 vg = c_light/(n1+(n1-n0)/std::log(E1/E0));
636 <<
637 // allow only for 'normal dispersion' -> d 278 // allow only for 'normal dispersion' -> dn/d(logE) > 0
638 if ((vg < 0) || (vg > c_light / n1)) { << 279 if(vg<0 || vg>c_light/n1)vg = c_light/n1;
639 vg = c_light / n1; << 280 groupvel->AddElement( E1, vg );
640 } << 281 }else{
641 groupvel->InsertValues(E1, vg); << 282 // only one entry in RINDEX -- weird!
>> 283 groupvel->AddElement( E0, c_light/n0 );
642 } 284 }
643 else // only one entry in RINDEX -- weird! << 285
644 { <<
645 groupvel->InsertValues(E0, c_light / n0); <<
646 } <<
647 <<
648 this->AddProperty("GROUPVEL", groupvel); <<
649 <<
650 return groupvel; 286 return groupvel;
651 } 287 }
652 288