1 | // |
---|
2 | // ******************************************************************** |
---|
3 | // * License and Disclaimer * |
---|
4 | // * * |
---|
5 | // * The Geant4 software is copyright of the Copyright Holders of * |
---|
6 | // * the Geant4 Collaboration. It is provided under the terms and * |
---|
7 | // * conditions of the Geant4 Software License, included in the file * |
---|
8 | // * LICENSE and available at http://cern.ch/geant4/license . These * |
---|
9 | // * include a list of copyright holders. * |
---|
10 | // * * |
---|
11 | // * Neither the authors of this software system, nor their employing * |
---|
12 | // * institutes,nor the agencies providing financial support for this * |
---|
13 | // * work make any representation or warranty, express or implied, * |
---|
14 | // * regarding this software system or assume any liability for its * |
---|
15 | // * use. Please see the license in the file LICENSE and URL above * |
---|
16 | // * for the full disclaimer and the limitation of liability. * |
---|
17 | // * * |
---|
18 | // * This code implementation is the result of the scientific and * |
---|
19 | // * technical work of the GEANT4 collaboration. * |
---|
20 | // * By using, copying, modifying or distributing the software (or * |
---|
21 | // * any work based on the software) you agree to acknowledge its * |
---|
22 | // * use in resulting scientific publications, and indicate your * |
---|
23 | // * acceptance of all terms of the Geant4 Software license. * |
---|
24 | // ******************************************************************** |
---|
25 | // |
---|
26 | //////////////////////////////////////////////////////////////////////////////// |
---|
27 | // |
---|
28 | #include "exrdmMaterial.hh" |
---|
29 | #include "exrdmMaterialData.hh" |
---|
30 | #include "exrdmMaterialMessenger.hh" |
---|
31 | |
---|
32 | #include "globals.hh" |
---|
33 | #include "G4UnitsTable.hh" |
---|
34 | #include "G4ios.hh" |
---|
35 | #include <vector> |
---|
36 | #include <iomanip> |
---|
37 | //////////////////////////////////////////////////////////////////////////////// |
---|
38 | // |
---|
39 | exrdmMaterial::exrdmMaterial () |
---|
40 | { |
---|
41 | Material.clear(); |
---|
42 | Element.clear(); |
---|
43 | Isotope.clear(); |
---|
44 | // some default materials vacuum (0), air (1) and aluminium (2) defined here |
---|
45 | // examples of vacuum |
---|
46 | // |
---|
47 | // G4double a,z; |
---|
48 | |
---|
49 | static G4bool bmat = false ; |
---|
50 | |
---|
51 | if (!bmat) { |
---|
52 | |
---|
53 | // vacuum |
---|
54 | G4double density = universe_mean_density; //from PhysicalConstants.h |
---|
55 | G4double pressure = 3.e-18*pascal; |
---|
56 | G4double temperature = 2.73*kelvin; |
---|
57 | AddMaterial("Vacuum", "H", density,"gas",temperature,pressure); |
---|
58 | // air |
---|
59 | density = 1.290*mg/cm3; |
---|
60 | AddMaterial("Air", "N0.78-O0.22", density, "gas"); |
---|
61 | |
---|
62 | // aluminium |
---|
63 | density=2.700*g/cm3 ; |
---|
64 | AddMaterial ("Aluminium", "Al", density,""); |
---|
65 | |
---|
66 | //silicon |
---|
67 | density=2.3290*g/cm3 ; |
---|
68 | AddMaterial ("Silicon", "Si", density,""); |
---|
69 | |
---|
70 | bmat = true; |
---|
71 | } |
---|
72 | // create commands for interactive definition of the geometry |
---|
73 | materialMessenger = new exrdmMaterialMessenger(this); |
---|
74 | } |
---|
75 | //////////////////////////////////////////////////////////////////////////////// |
---|
76 | // |
---|
77 | exrdmMaterial::~exrdmMaterial () |
---|
78 | { |
---|
79 | delete materialMessenger; |
---|
80 | } |
---|
81 | //////////////////////////////////////////////////////////////////////////////// |
---|
82 | // |
---|
83 | void exrdmMaterial::AddMaterial (G4String name, G4String formula, G4double density, |
---|
84 | G4String state, G4double tem, G4double pres) |
---|
85 | { |
---|
86 | G4int isotope, Z; |
---|
87 | size_t i; |
---|
88 | for (i = 0; i<Material.size(); i++) { |
---|
89 | if (Material[i]->GetName() == name) { |
---|
90 | G4cerr <<" AddMaterial : material " <<name |
---|
91 | <<" already exists." <<G4endl; |
---|
92 | G4cerr <<"--> Command rejected." <<G4endl; |
---|
93 | return; |
---|
94 | } |
---|
95 | } |
---|
96 | |
---|
97 | char *tokenPtr1 = NULL; |
---|
98 | char *sname = NULL; |
---|
99 | G4String s, s1("0123456789"); |
---|
100 | G4String element, isotopename; |
---|
101 | G4int ncomponents, natoms; |
---|
102 | G4double fatoms = 0.; |
---|
103 | size_t ls, id=0, ll, lr; |
---|
104 | ncomponents = 0; |
---|
105 | |
---|
106 | sname = new char[strlen(formula)+1]; |
---|
107 | strcpy(sname,formula); |
---|
108 | tokenPtr1 = strtok(sname,"-"); |
---|
109 | |
---|
110 | while (tokenPtr1 != NULL) { |
---|
111 | ncomponents++; |
---|
112 | tokenPtr1 = strtok( NULL, "-"); |
---|
113 | } |
---|
114 | delete[] sname; |
---|
115 | |
---|
116 | G4Material* aMaterial = 0; |
---|
117 | // G4cout << name <<" "<< formula << " " << density/(g/cm3) << " " << tem <<" " <<pres << G4endl; |
---|
118 | |
---|
119 | if (state == "") { |
---|
120 | aMaterial = new G4Material(name, density, ncomponents); |
---|
121 | } else if (state == "solid" && tem > 0.) { |
---|
122 | aMaterial = new G4Material(name, density, ncomponents, |
---|
123 | kStateSolid, tem ); |
---|
124 | } else if (state == "gas" && pres > 0.) { |
---|
125 | aMaterial = new G4Material(name, density, ncomponents, |
---|
126 | kStateGas, tem, pres ); |
---|
127 | } |
---|
128 | if (aMaterial == 0) { |
---|
129 | G4cerr <<" AddMaterial : Name " <<name <<"." <<G4endl; |
---|
130 | G4cerr <<"--> Command failed." <<G4endl; |
---|
131 | return; |
---|
132 | } |
---|
133 | |
---|
134 | sname=new char[strlen(formula)+1]; |
---|
135 | strcpy(sname,formula); |
---|
136 | tokenPtr1 = strtok(sname,"-"); |
---|
137 | |
---|
138 | while (tokenPtr1 != NULL) { |
---|
139 | isotope = 0; |
---|
140 | // G4cout << tokenPtr1 << G4endl; |
---|
141 | s = G4String(tokenPtr1); |
---|
142 | ls = s.length(); |
---|
143 | ll = s.find("("); |
---|
144 | lr = s.find(")"); |
---|
145 | if (ll == lr) { |
---|
146 | id = s.find_first_of(s1); |
---|
147 | element = s.substr(0,id); |
---|
148 | |
---|
149 | if (element.length() == 1) element.insert(0," "); |
---|
150 | for (i = 0; i<110; i++) { |
---|
151 | if (element == ELU[i]) break; |
---|
152 | } |
---|
153 | if (i == 110) { |
---|
154 | for (i = 0; i<110; i++) { |
---|
155 | if (element == ELL[i]) break; |
---|
156 | } |
---|
157 | if (i == 110) { |
---|
158 | for (i = 0; i<110; i++) { |
---|
159 | if (element == EUU[i]) break; |
---|
160 | } |
---|
161 | } |
---|
162 | } |
---|
163 | |
---|
164 | if (i == 110) { |
---|
165 | G4cerr <<"AddMaterial : Invalid element in material formula." |
---|
166 | <<element <<G4endl; |
---|
167 | G4cerr <<"--> Command rejected." <<G4endl; |
---|
168 | // delete aMaterial; |
---|
169 | // Material[NbMat] = NULL; |
---|
170 | return; |
---|
171 | } |
---|
172 | |
---|
173 | Z = i+1; |
---|
174 | element = ELU[i]; |
---|
175 | if (id == std::string::npos) { |
---|
176 | natoms = 1; |
---|
177 | } else { |
---|
178 | natoms = atoi((s.substr(id,ls-id)).c_str()); |
---|
179 | } |
---|
180 | if (natoms < 1) fatoms = atof((s.substr(id,ls-id)).c_str()); |
---|
181 | // G4cout << " Elements = " << element << G4endl; |
---|
182 | //G4cout << " Nb of atoms = " << natoms << G4endl; |
---|
183 | } else { |
---|
184 | element = s.substr(0,ll); |
---|
185 | isotope = atoi((s.substr(ll+1,lr-ll)).c_str()); |
---|
186 | if (element.length() == 1) element.insert(0," "); |
---|
187 | for (i = 0; i<110; i++) { |
---|
188 | if (element == ELU[i]) break; |
---|
189 | } |
---|
190 | if (i == 110) { |
---|
191 | for (i = 0; i<110; i++) { |
---|
192 | if (element == ELL[i]) break; |
---|
193 | } |
---|
194 | if (i == 110) { |
---|
195 | for (i = 0; i<110; i++) { |
---|
196 | if (element == EUU[i]) break; |
---|
197 | } |
---|
198 | } |
---|
199 | } |
---|
200 | if (i == 110) { |
---|
201 | G4cerr <<"AddMaterial : Invalid element in material formula." |
---|
202 | <<element <<G4endl; |
---|
203 | G4cerr <<"--> Command rejected." <<G4endl; |
---|
204 | // delete aMaterial; |
---|
205 | // Material[NbMat] = NULL; |
---|
206 | return; |
---|
207 | } |
---|
208 | |
---|
209 | Z = i+1; |
---|
210 | element = ELU[i]; |
---|
211 | isotopename = element+s.substr(ll+1,lr-ll-1); |
---|
212 | if (lr == std::string::npos ) { |
---|
213 | natoms = 1; |
---|
214 | } else { |
---|
215 | natoms = atoi((s.substr(lr+1,ls-lr)).c_str()); |
---|
216 | } |
---|
217 | if (natoms < 1) fatoms = atof((s.substr(id,ls-id)).c_str()); |
---|
218 | if (fatoms == 0.) natoms = 1; |
---|
219 | // |
---|
220 | // G4cout << " Elements = " << element << G4endl; |
---|
221 | // G4cout << " Isotope Nb = " << isotope << G4endl; |
---|
222 | // G4cout << " Nb of atoms = " << natoms << G4endl; |
---|
223 | } |
---|
224 | if (isotope != 0) { |
---|
225 | if (G4Isotope::GetIsotope(isotopename) == NULL) { |
---|
226 | // G4Isotope* aIsotope = new G4Isotope(isotopename, Z, isotope, A[Z-1]*g/mole); |
---|
227 | G4Isotope* aIsotope = new G4Isotope(isotopename, Z, isotope, isotope*g/mole); |
---|
228 | G4Element* aElement = new G4Element(isotopename, element, 1); |
---|
229 | aElement->AddIsotope(aIsotope, 100.*perCent); |
---|
230 | Isotope.push_back(aIsotope); |
---|
231 | if (natoms>0) { |
---|
232 | aMaterial->AddElement(aElement, natoms); |
---|
233 | } else { |
---|
234 | aMaterial->AddElement(aElement, fatoms); |
---|
235 | } |
---|
236 | Element.push_back(aElement); |
---|
237 | } else { |
---|
238 | if (natoms>0) { |
---|
239 | aMaterial->AddElement( G4Element::GetElement(isotopename) , natoms); |
---|
240 | } else { |
---|
241 | aMaterial->AddElement( G4Element::GetElement(isotopename) , fatoms); |
---|
242 | } |
---|
243 | } |
---|
244 | } else { |
---|
245 | if ( G4Element::GetElement(element) == NULL) { |
---|
246 | G4Element* aElement = new G4Element(element, element, Z, A[Z-1]*g/mole); |
---|
247 | if (natoms>0) { |
---|
248 | aMaterial->AddElement(aElement, natoms); |
---|
249 | } else { |
---|
250 | aMaterial->AddElement(aElement, fatoms); |
---|
251 | } |
---|
252 | Element.push_back(aElement); |
---|
253 | } else { |
---|
254 | if (natoms>0) { |
---|
255 | aMaterial->AddElement( G4Element::GetElement(element) , natoms); |
---|
256 | } else { |
---|
257 | aMaterial->AddElement( G4Element::GetElement(element) , fatoms); |
---|
258 | } |
---|
259 | } |
---|
260 | } |
---|
261 | tokenPtr1 = strtok( NULL, "-"); |
---|
262 | // s.empty(); |
---|
263 | //element.erase(); |
---|
264 | // |
---|
265 | } |
---|
266 | |
---|
267 | delete[] sname; |
---|
268 | Material.push_back(aMaterial); |
---|
269 | G4cout <<" Material:" <<name <<" with formula: " <<formula <<" added! " |
---|
270 | <<G4endl; |
---|
271 | G4cout <<" Nb of Material = " <<Material.size() <<G4endl; |
---|
272 | G4cout <<" Nb of Isotope = " <<Isotope.size() <<G4endl; |
---|
273 | G4cout <<" Nb of Element = " <<Element.size() <<G4endl; |
---|
274 | } |
---|
275 | //////////////////////////////////////////////////////////////////////////////// |
---|
276 | // |
---|
277 | void exrdmMaterial::DeleteMaterial (G4int j) |
---|
278 | { |
---|
279 | size_t i(j-1); |
---|
280 | if (i > Material.size()) { |
---|
281 | G4cerr <<"DeleteMaterial : Invalid material index " <<j <<"." <<G4endl; |
---|
282 | G4cerr <<"--> Command rejected." <<G4endl; |
---|
283 | } else { |
---|
284 | G4cerr <<"It seems there is no mechanism in G4 for deleting a material yet!" |
---|
285 | <<G4endl; |
---|
286 | G4cerr <<"--> Command rejected." <<G4endl; |
---|
287 | } |
---|
288 | } |
---|
289 | //////////////////////////////////////////////////////////////////////////////// |
---|
290 | // |
---|
291 | void exrdmMaterial::DeleteMaterial (G4String ) |
---|
292 | { |
---|
293 | G4cerr <<"It seems there is no mechanism in G4 for deleting a material yet!" |
---|
294 | <<G4endl; |
---|
295 | G4cerr <<"--> Command rejected." <<G4endl; |
---|
296 | } |
---|
297 | //////////////////////////////////////////////////////////////////////////////// |
---|
298 | // |
---|
299 | G4int exrdmMaterial::GetMaterialIndex (G4String name) |
---|
300 | { |
---|
301 | size_t i ; |
---|
302 | for (i = 0; i < Material.size(); i++) { |
---|
303 | if (Material[i]->GetName() == name) break; |
---|
304 | } |
---|
305 | G4int k = G4int(i); |
---|
306 | if (i == Material.size()) k = -1; |
---|
307 | return k; |
---|
308 | } |
---|
309 | //////////////////////////////////////////////////////////////////////////////// |
---|
310 | // |
---|
311 | void exrdmMaterial::ListMaterial () |
---|
312 | { |
---|
313 | G4cout <<" There are" <<std::setw(3) <<Material.size() |
---|
314 | <<" materials defined." <<G4endl; |
---|
315 | for (size_t i = 0; i< Material.size(); i++) |
---|
316 | G4cout <<" Material Index " <<std::setw(3) <<i+1 <<" " |
---|
317 | <<std::setw(14) <<Material[i]->GetName() |
---|
318 | <<" density: " <<std::setw(6) <<std::setprecision(3) |
---|
319 | <<G4BestUnit(Material[i]->GetDensity(),"Volumic Mass") <<G4endl; |
---|
320 | G4cout <<G4endl; |
---|
321 | |
---|
322 | } |
---|
323 | //////////////////////////////////////////////////////////////////////////////// |
---|