1 | |
2 | |
3 | |
4 | |
5 | |
6 | |
7 | |
8 | #define APP_NAME"HddsG4Builder" "HddsG4Builder" |
9 | |
10 | #include <HddsG4Builder.hh> |
11 | #include <G4Box.hh> |
12 | #include <G4Trd.hh> |
13 | #include <G4Trap.hh> |
14 | #include <G4Para.hh> |
15 | #include <G4Tubs.hh> |
16 | #include <G4Cons.hh> |
17 | #include <G4Hype.hh> |
18 | #include <G4Torus.hh> |
19 | #include <G4Orb.hh> |
20 | #include <G4Sphere.hh> |
21 | #include <G4Polycone.hh> |
22 | #include <G4Polyhedra.hh> |
23 | #include <G4EllipticalTube.hh> |
24 | #include <G4SystemOfUnits.hh> |
25 | #include <G4Mag_UsualEqRhs.hh> |
26 | #include <G4MagIntegratorStepper.hh> |
27 | #include <G4ExactHelixStepper.hh> |
28 | #include <G4HelixMixedStepper.hh> |
29 | #include <G4ClassicalRK4.hh> |
30 | #include <G4ChordFinder.hh> |
31 | #include <G4PVPlacement.hh> |
32 | #include <G4PVDivision.hh> |
33 | #include <G4OpticalSurface.hh> |
34 | #include <G4LogicalSkinSurface.hh> |
35 | #include <G4VisAttributes.hh> |
36 | |
37 | #include <xercesc/util/PlatformUtils.hpp> |
38 | #include <xercesc/util/XMLString.hpp> |
39 | #include <xercesc/util/XMLStringTokenizer.hpp> |
40 | #include <xercesc/sax/SAXParseException.hpp> |
41 | #include <xercesc/parsers/XercesDOMParser.hpp> |
42 | #include <xercesc/framework/LocalFileFormatTarget.hpp> |
43 | #include <xercesc/dom/DOM.hpp> |
44 | #include <xercesc/util/XercesDefs.hpp> |
45 | #include <xercesc/sax/ErrorHandler.hpp> |
46 | |
47 | using namespace xercesc; |
48 | |
49 | #include <assert.h> |
50 | #include <stdlib.h> |
51 | #include <ctype.h> |
52 | #include <stdio.h> |
53 | #include <math.h> |
54 | |
55 | #include <iostream> |
56 | #include <fstream> |
57 | #include <sstream> |
58 | #include <string> |
59 | #include <iomanip> |
60 | #include <vector> |
61 | #include <list> |
62 | #include <map> |
63 | |
64 | #define X(str)XString(str).unicode_str() XString(str).unicode_str() |
65 | #define S(str)str.c_str() str.c_str() |
66 | |
67 | #ifdef LINUX_CPUTIME_PROFILING |
68 | extern CPUtimer timer; |
69 | #endif |
70 | |
71 | HddsG4Builder::HddsG4Builder() : fWorldVolume(0) { } |
72 | |
73 | HddsG4Builder::HddsG4Builder(const HddsG4Builder &src) |
74 | { |
75 | fWorldVolume = src.fWorldVolume; |
76 | fElements = src.fElements; |
77 | fMaterials = src.fMaterials; |
78 | fMagneticRegions = src.fMagneticRegions; |
79 | fLogicalVolumes = src.fLogicalVolumes; |
80 | fPhysicalVolumes = src.fPhysicalVolumes; |
81 | fSensitiveVolumes = src.fSensitiveVolumes; |
82 | fRotations = src.fRotations; |
83 | fCurrentMother = src.fCurrentMother; |
84 | fCurrentPlacement = src.fCurrentPlacement; |
85 | fCurrentDivision = src.fCurrentDivision; |
86 | fCurrentPhiCenter = src.fCurrentPhiCenter; |
87 | } |
88 | |
89 | HddsG4Builder::~HddsG4Builder() { } |
90 | |
91 | int HddsG4Builder::createMaterial(DOMElement* el) |
92 | { |
93 | #ifdef LINUX_CPUTIME_PROFILING |
94 | std::stringstream timestr; |
95 | timestr << "createMaterial: " << timer.getUserTime() << ", " |
96 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
97 | timer.resetClocks(); |
98 | #endif |
99 | int imate = CodeWriter::createMaterial(el); |
100 | |
101 | if (fSubst.fBrewList.size() == 0) |
102 | { |
103 | XString matS = fSubst.getName(); |
104 | XString symS = fSubst.getSymbol(); |
105 | double A = fSubst.getAtomicWeight(); |
106 | double Z = fSubst.getAtomicNumber(); |
107 | double dens = fSubst.getDensity(); |
108 | fElements[imate] = new G4Element(matS,symS,Z,A*g/mole); |
109 | fMaterials[imate] = new G4Material(matS,Z,A*g/mole,dens*g/cm3); |
110 | } |
111 | else |
112 | { |
113 | XString matS = fSubst.getName(); |
114 | double dens = fSubst.getDensity(); |
115 | int ncomp = fSubst.fBrewList.size(); |
116 | fMaterials[imate] = new G4Material(matS,dens*g/cm3,ncomp); |
117 | std::list<Substance::Brew>::iterator iter; |
118 | for (iter = fSubst.fBrewList.begin(); |
119 | iter != fSubst.fBrewList.end(); iter++) |
120 | { |
121 | int subimate = iter->sub->fUniqueID; |
122 | std::map<int,G4Element*>::iterator elemit = fElements.find(subimate); |
123 | if (iter->natoms) |
124 | { |
125 | fMaterials[imate]->AddElement(elemit->second,iter->natoms); |
126 | } |
127 | else if (elemit != fElements.end()) |
128 | { |
129 | fMaterials[imate]->AddElement(elemit->second,iter->wfact); |
130 | } |
131 | else |
132 | { |
133 | fMaterials[imate]->AddMaterial(fMaterials[subimate],iter->wfact); |
134 | } |
135 | } |
136 | } |
137 | |
138 | DOMNodeList* propList = el->getElementsByTagName(X("optical_properties")XString("optical_properties").unicode_str()); |
139 | if (propList->getLength() > 0) |
140 | { |
141 | DOMElement* propEl = (DOMElement*)propList->item(0); |
142 | DOMNodeList* specL = propEl->getElementsByTagName(X("specify")XString("specify").unicode_str()); |
143 | int len = specL->getLength(); |
144 | std::vector<double> Ephot; |
145 | std::vector<double> abslen; |
146 | std::vector<double> effic; |
147 | std::vector<double> rindex; |
148 | std::vector<double> smooth; |
149 | std::vector<double> reflect; |
150 | for (int i=0; i < len; ++i) |
151 | { |
152 | DOMElement* specEl = (DOMElement*)specL->item(i); |
153 | XString valS; |
154 | valS = specEl->getAttribute(X("E")XString("E").unicode_str()); |
155 | Ephot.push_back(atof(S(valS)valS.c_str())*eV); |
156 | valS = specEl->getAttribute(X("refindex")XString("refindex").unicode_str()); |
157 | rindex.push_back(atof(S(valS)valS.c_str())); |
158 | valS = specEl->getAttribute(X("abslen")XString("abslen").unicode_str()); |
159 | abslen.push_back(atof(S(valS)valS.c_str())*cm); |
160 | valS = specEl->getAttribute(X("smooth")XString("smooth").unicode_str()); |
161 | smooth.push_back(atof(S(valS)valS.c_str())); |
162 | valS = specEl->getAttribute(X("reflect")XString("reflect").unicode_str()); |
163 | reflect.push_back(atof(S(valS)valS.c_str())); |
164 | valS = specEl->getAttribute(X("effic")XString("effic").unicode_str()); |
165 | effic.push_back(atof(S(valS)valS.c_str())); |
166 | } |
167 | G4MaterialPropertiesTable *mpt = new G4MaterialPropertiesTable(); |
168 | if (rindex[0] > 0) |
169 | { |
170 | mpt->AddProperty("RINDEX",&Ephot[0],&rindex[0],len); |
171 | } |
172 | if (abslen[0] > 0) |
173 | { |
174 | mpt->AddProperty("ABSLENGTH",&Ephot[0],&abslen[0],len); |
175 | } |
176 | if (smooth[0] > 0) |
177 | { |
178 | mpt->AddProperty("POLISH",&Ephot[0],&smooth[0],len); |
179 | } |
180 | if (reflect[0] > 0) |
181 | { |
182 | mpt->AddProperty("REFLECTIVITY",&Ephot[0],&reflect[0],len); |
183 | } |
184 | if (effic[0] > 0) |
185 | { |
186 | mpt->AddProperty("EFFICIENCY",&Ephot[0],&effic[0],len); |
187 | } |
188 | fMaterials[imate]->SetMaterialPropertiesTable(mpt); |
189 | } |
190 | #ifdef LINUX_CPUTIME_PROFILING |
191 | timestr << " ( " << timer.getUserDelta() << " ) "; |
192 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
193 | #endif |
194 | return imate; |
195 | } |
196 | |
197 | int HddsG4Builder::createSolid(DOMElement* el, Refsys& ref) |
198 | { |
199 | #ifdef LINUX_CPUTIME_PROFILING |
200 | std::stringstream timestr; |
201 | timestr << "createSolid: " << timer.getUserTime() << ", " |
202 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
203 | timer.resetClocks(); |
204 | #endif |
205 | |
206 | int ivolu = CodeWriter::createSolid(el,ref); |
207 | XString nameS(el->getAttribute(X("name")XString("name").unicode_str())); |
208 | int imate = fSubst.fUniqueID; |
209 | |
210 | Units unit; |
211 | unit.getConversions(el); |
212 | |
213 | G4VSolid *solid; |
214 | XString shapeS(el->getTagName()); |
215 | if (shapeS == "box") |
216 | { |
217 | double xl, yl, zl; |
218 | XString xyzS(el->getAttribute(X("X_Y_Z")XString("X_Y_Z").unicode_str())); |
219 | std::stringstream listr(xyzS); |
220 | listr >> xl >> yl >> zl; |
221 | solid = new G4Box(S(nameS)nameS.c_str(), |
222 | xl/2 *cm/unit.cm, |
223 | yl/2 *cm/unit.cm, |
224 | zl/2 *cm/unit.cm |
225 | ); |
226 | } |
227 | else if (shapeS == "tubs") |
228 | { |
229 | double ri, ro, zl, phi0, dphi; |
230 | XString riozS(el->getAttribute(X("Rio_Z")XString("Rio_Z").unicode_str())); |
231 | std::stringstream listr(riozS); |
232 | listr >> ri >> ro >> zl; |
233 | XString profS(el->getAttribute(X("profile")XString("profile").unicode_str())); |
234 | listr.clear(), listr.str(profS); |
235 | listr >> phi0 >> dphi; |
236 | solid = new G4Tubs(S(nameS)nameS.c_str(), |
237 | ri * cm/unit.cm, |
238 | ro * cm/unit.cm, |
239 | zl/2 * cm/unit.cm, |
240 | phi0 * deg/unit.deg, |
241 | dphi * deg/unit.deg |
242 | ); |
243 | } |
244 | else if (shapeS == "eltu") |
245 | { |
246 | |
247 | double rx, ry, zl; |
248 | XString rxyzS(el->getAttribute(X("Rxy_Z")XString("Rxy_Z").unicode_str())); |
249 | std::stringstream listr(rxyzS); |
250 | listr >> rx >> ry >> zl; |
251 | solid = new G4EllipticalTube(S(nameS)nameS.c_str(), |
252 | rx * cm/unit.cm, |
253 | ry * cm/unit.cm, |
254 | zl/2 * cm/unit.cm |
255 | ); |
256 | } |
257 | else if (shapeS == "trd") |
258 | { |
259 | double xm, ym, xp, yp, zl; |
260 | XString xyzS(el->getAttribute(X("Xmp_Ymp_Z")XString("Xmp_Ymp_Z").unicode_str())); |
261 | std::stringstream listr(xyzS); |
262 | listr >> xm >> xp >> ym >> yp >> zl; |
263 | xm += 1e-100; |
264 | xp += 1e-100; |
265 | ym += 1e-100; |
266 | yp += 1e-100; |
267 | zl += 1e-100; |
268 | double alph_xz, alph_yz; |
269 | XString incS(el->getAttribute(X("inclination")XString("inclination").unicode_str())); |
270 | listr.clear(), listr.str(incS); |
271 | listr >> alph_xz >> alph_yz; |
272 | double x = tan(alph_xz/unit.rad); |
273 | double y = tan(alph_yz/unit.rad); |
274 | double r = sqrt(x*x + y*y); |
275 | solid = new G4Trap(S(nameS)nameS.c_str(), |
276 | zl/2 * cm/unit.cm, |
277 | atan2(r,1) * rad, |
278 | atan2(y,x) * rad, |
279 | ym/2 * cm/unit.cm, |
280 | xm/2 * cm/unit.cm, |
281 | xm/2 * cm/unit.cm, |
282 | 0, |
283 | yp/2 * cm/unit.cm, |
284 | xp/2 * cm/unit.cm, |
285 | xp/2 * cm/unit.cm, |
286 | 0 |
287 | ); |
288 | } |
289 | else if (shapeS == "pcon") |
290 | { |
291 | double phi0, dphi; |
292 | XString profS(el->getAttribute(X("profile")XString("profile").unicode_str())); |
293 | std::stringstream listr(profS); |
294 | listr >> phi0 >> dphi; |
295 | DOMNodeList* planeList = el->getElementsByTagName(X("polyplane")XString("polyplane").unicode_str()); |
296 | int nplanes = planeList->getLength(); |
297 | std::vector<double> zPlane; |
298 | std::vector<double> rInner; |
299 | std::vector<double> rOuter; |
300 | double zlast = -1e30; |
301 | double zeps = 0; |
302 | for (unsigned int p = 0; p < planeList->getLength(); p++) |
303 | { |
304 | double ri, ro, zl; |
305 | DOMNode* node = planeList->item(p); |
306 | DOMElement* elem = (DOMElement*) node; |
307 | XString riozS(elem->getAttribute(X("Rio_Z")XString("Rio_Z").unicode_str())); |
308 | std::stringstream listr1(riozS); |
309 | listr1 >> ri >> ro >> zl; |
310 | if (zl < zlast) |
311 | { |
312 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: Please re-order the polyplanes" |
313 | << " of volume " << S(nameS)nameS.c_str() << " so that the z-values" |
314 | << " are non-decreasing" |
315 | << G4endlstd::endl; |
316 | exit(1); |
317 | } |
318 | zl = (zl < zlast + zeps)? zlast + zeps : zl; |
319 | zlast = zl; |
320 | zPlane.push_back(zl * cm/unit.cm); |
321 | rInner.push_back(ri * cm/unit.cm); |
322 | rOuter.push_back(ro * cm/unit.cm); |
323 | } |
324 | solid = new G4Polycone(S(nameS)nameS.c_str(), |
325 | phi0 * deg/unit.deg, |
326 | dphi * deg/unit.deg, |
327 | nplanes, |
328 | &zPlane[0], |
329 | &rInner[0], |
330 | &rOuter[0] |
331 | ); |
332 | } |
333 | else if (shapeS == "pgon") |
334 | { |
335 | int segments; |
336 | XString segS(el->getAttribute(X("segments")XString("segments").unicode_str())); |
337 | segments = atoi(S(segS)segS.c_str()); |
338 | double phi0, dphi; |
339 | XString profS(el->getAttribute(X("profile")XString("profile").unicode_str())); |
340 | std::stringstream listr(profS); |
341 | listr >> phi0 >> dphi; |
342 | DOMNodeList* planeList = el->getElementsByTagName(X("polyplane")XString("polyplane").unicode_str()); |
343 | int nplanes = planeList->getLength(); |
344 | std::vector<double> zPlane; |
345 | std::vector<double> rInner; |
346 | std::vector<double> rOuter; |
347 | double zlast = -1e30; |
348 | double zeps = 0; |
349 | for (unsigned int p = 0; p < planeList->getLength(); p++) |
350 | { |
351 | double ri, ro, zl; |
352 | DOMNode* node = planeList->item(p); |
353 | DOMElement* elem = (DOMElement*) node; |
354 | XString riozS(elem->getAttribute(X("Rio_Z")XString("Rio_Z").unicode_str())); |
355 | std::stringstream listr1(riozS); |
356 | listr1 >> ri >> ro >> zl; |
357 | if (zl < zlast) |
358 | { |
359 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: Please re-order the polyplanes" |
360 | << " of volume " << S(nameS)nameS.c_str() << " so that the z-values" |
361 | << " are non-decreasing." |
362 | << G4endlstd::endl; |
363 | exit(1); |
364 | } |
365 | zl = (zl < zlast + zeps)? zlast + zeps : zl; |
366 | zlast = zl; |
367 | zPlane.push_back(zl * cm/unit.cm); |
368 | rInner.push_back(ri * cm/unit.cm); |
369 | rOuter.push_back(ro * cm/unit.cm); |
370 | } |
371 | solid = new G4Polyhedra(S(nameS)nameS.c_str(), |
372 | phi0 * deg/unit.deg, |
373 | dphi * deg/unit.deg, |
374 | segments, |
375 | nplanes, |
376 | &zPlane[0], |
377 | &rInner[0], |
378 | &rOuter[0] |
379 | ); |
380 | } |
381 | else if (shapeS == "cons") |
382 | { |
383 | double rim, rip, rom, rop, zl; |
384 | XString riozS(el->getAttribute(X("Rio1_Rio2_Z")XString("Rio1_Rio2_Z").unicode_str())); |
385 | std::stringstream listr(riozS); |
386 | listr >> rim >> rom >> rip >> rop >> zl; |
387 | double phi0, dphi; |
388 | XString profS(el->getAttribute(X("profile")XString("profile").unicode_str())); |
389 | listr.clear(), listr.str(profS); |
390 | listr >> phi0 >> dphi; |
391 | solid = new G4Cons(S(nameS)nameS.c_str(), |
392 | zl/2 * cm/unit.cm, |
393 | rim * cm/unit.cm, |
394 | rom * cm/unit.cm, |
395 | rip * cm/unit.cm, |
396 | rop * cm/unit.cm, |
397 | phi0 * deg/unit.deg, |
398 | dphi * deg/unit.deg |
399 | ); |
400 | } |
401 | else if (shapeS == "sphere") |
402 | { |
403 | double ri, ro; |
404 | XString rioS(el->getAttribute(X("Rio")XString("Rio").unicode_str())); |
405 | std::stringstream listr(rioS); |
406 | listr >> ri >> ro; |
407 | double theta0, theta1; |
408 | XString polarS(el->getAttribute(X("polar_bounds")XString("polar_bounds").unicode_str())); |
409 | listr.clear(), listr.str(polarS); |
410 | listr >> theta0 >> theta1; |
411 | double phi0, dphi; |
412 | XString profS(el->getAttribute(X("profile")XString("profile").unicode_str())); |
413 | listr.clear(), listr.str(profS); |
414 | listr >> phi0 >> dphi; |
415 | solid = new G4Sphere(S(nameS)nameS.c_str(), |
416 | ri * cm/unit.cm, |
417 | ro * cm/unit.cm, |
418 | phi0 * deg/unit.deg, |
419 | dphi * deg/unit.deg, |
420 | theta0 * deg/unit.deg, |
421 | theta1 * deg/unit.deg |
422 | ); |
423 | } |
424 | else |
425 | { |
426 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: volume " << S(nameS)nameS.c_str() |
427 | << " should be one of the valid shapes, not " << S(shapeS)shapeS.c_str() |
428 | << G4endlstd::endl; |
429 | exit(1); |
430 | } |
431 | |
432 | vpair_t newvol(ivolu,0); |
433 | G4FieldManager* fieldmgr = fFieldManagers[ref.fRegionID]; |
434 | fLogicalVolumes[newvol] = new G4LogicalVolume(solid,fMaterials[imate], |
435 | S(nameS)nameS.c_str(),fieldmgr); |
436 | G4Material* mate = fMaterials[imate]; |
437 | double dens = mate->GetDensity(); |
438 | double dmod = (int)(dens*97.345/(g/cm3)) % 20; |
439 | double red = 1-3/(3+dens /(g/cm3)); |
440 | double green = 1-5/(5+dmod); |
441 | double blue = 1-green; |
442 | G4VisAttributes* attr = new G4VisAttributes(G4Colour(red,green,blue)); |
443 | fLogicalVolumes[newvol]->SetVisAttributes(attr); |
444 | |
445 | XString sensiS(el->getAttribute(X("sensitive")XString("sensitive").unicode_str())); |
446 | if (sensiS == "true") |
447 | { |
448 | fSensitiveVolumes[ivolu] = fLogicalVolumes[newvol]; |
449 | } |
450 | |
451 | G4MaterialPropertiesTable *mpt = fMaterials[imate]->GetMaterialPropertiesTable(); |
452 | if (mpt) |
453 | { |
454 | G4MaterialPropertyVector *poli_vector = mpt->GetProperty("POLISH"); |
455 | G4MaterialPropertyVector *refl_vector = mpt->GetProperty("REFLECTIVITY"); |
456 | G4MaterialPropertyVector *refi_vector = mpt->GetProperty("RINDEX"); |
457 | G4MaterialPropertyVector *absl_vector = mpt->GetProperty("ABSLENGTH"); |
458 | G4MaterialPropertyVector *effi_vector = mpt->GetProperty("EFFICIENCY"); |
459 | if (refi_vector == 0) { |
460 | if (refl_vector == 0) { |
461 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: material " << mate->GetName() |
462 | << " needs to have either refindex defined (dielectrics)" |
463 | << " or else reflect (metals) but it has neither!" |
464 | << G4endlstd::endl; |
465 | exit(1); |
466 | } |
467 | if (absl_vector != 0) { |
468 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: material " << mate->GetName() |
469 | << " has abslen defined but no refindex; defining abslen" |
470 | << " requires that refindex also be specified!" |
471 | << G4endlstd::endl; |
472 | exit(1); |
473 | } |
474 | if (poli_vector != 0 || effi_vector != 0) { |
475 | G4OpticalSurface *surface = new G4OpticalSurface(S(nameS)nameS.c_str()); |
476 | surface->SetType(dielectric_metal); |
477 | surface->SetModel(glisur); |
478 | surface->SetMaterialPropertiesTable(mpt); |
479 | if (poli_vector != 0) { |
480 | double polish = poli_vector->GetMaxValue(); |
481 | surface->SetPolish(polish); |
482 | surface->SetFinish(ground); |
483 | } |
484 | else { |
485 | surface->SetFinish(polished); |
486 | } |
487 | new G4LogicalSkinSurface(S(nameS)nameS.c_str(),fLogicalVolumes[newvol],surface); |
488 | } |
489 | } |
490 | else if (refl_vector == 0) { |
491 | if (poli_vector != 0) { |
492 | G4OpticalSurface *surface = new G4OpticalSurface(S(nameS)nameS.c_str()); |
493 | surface->SetType(dielectric_dielectric); |
494 | surface->SetModel(glisur); |
495 | surface->SetMaterialPropertiesTable(mpt); |
496 | if (poli_vector != 0) { |
497 | double polish = poli_vector->GetMaxValue(); |
498 | surface->SetPolish(polish); |
499 | surface->SetFinish(ground); |
500 | } |
501 | else { |
502 | surface->SetFinish(polished); |
503 | } |
504 | new G4LogicalSkinSurface(S(nameS)nameS.c_str(),fLogicalVolumes[newvol],surface); |
505 | } |
506 | } |
507 | else { |
508 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: material " << mate->GetName() |
509 | << " has optical properties, but neither refindex nor reflect" |
510 | << " data are provided, you need either one or the other!" |
511 | << G4endlstd::endl; |
512 | exit(1); |
513 | } |
514 | } |
515 | |
516 | #ifdef LINUX_CPUTIME_PROFILING |
517 | timestr << " ( " << timer.getUserDelta() << " ) "; |
518 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
519 | #endif |
520 | return ivolu; |
521 | } |
522 | |
523 | int HddsG4Builder::createRotation(Refsys& ref) |
524 | { |
525 | #ifdef LINUX_CPUTIME_PROFILING |
526 | std::stringstream timestr; |
527 | timestr << "createRotation: " << timer.getUserTime() << ", " |
528 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
529 | timer.resetClocks(); |
530 | #endif |
531 | int irot = CodeWriter::createRotation(ref); |
532 | |
533 | if (irot > 0) |
534 | { |
535 | std::vector<double> omega = ref.getRotation(); |
536 | fRotations[irot] = new G4RotationMatrix(); |
537 | fRotations[irot]->rotateX(omega[0]); |
538 | fRotations[irot]->rotateY(omega[1]); |
539 | fRotations[irot]->rotateZ(omega[2]); |
540 | fRotations[irot]->invert(); |
541 | } |
542 | else if (fRotations.count(0) == 0) |
543 | { |
544 | fRotations[0] = new G4RotationMatrix(); |
545 | } |
546 | |
547 | #ifdef LINUX_CPUTIME_PROFILING |
548 | timestr << " ( " << timer.getUserDelta() << " ) "; |
549 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
550 | #endif |
551 | return irot; |
552 | } |
553 | |
554 | int HddsG4Builder::createVolume(DOMElement* el, Refsys& ref) |
555 | { |
556 | #ifdef LINUX_CPUTIME_PROFILING |
557 | std::stringstream timestr; |
558 | timestr << "createVolume: " << timer.getUserTime() << ", " |
559 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
560 | timer.resetClocks(); |
561 | #endif |
562 | int icopy = CodeWriter::createVolume(el,ref); |
563 | |
564 | if (fPending) |
565 | { |
566 | XString myvoluS(el->getAttribute(X("HDDSvolu")XString("HDDSvolu").unicode_str())); |
567 | XString motherS(fRef.fMother->getAttribute(X("HDDSvolu")XString("HDDSvolu").unicode_str())); |
568 | int myvoluI = atoi(S(myvoluS)myvoluS.c_str()); |
569 | int motherI = atoi(S(motherS)motherS.c_str()); |
570 | G4ThreeVector origin(fRef.fOrigin[0]*cm, |
571 | fRef.fOrigin[1]*cm, |
572 | fRef.fOrigin[2]*cm); |
573 | int irot = fRef.fRotation; |
574 | |
575 | |
576 | |
577 | |
578 | G4RotationMatrix *rot = fRotations[irot]; |
579 | if (fCurrentPhiCenter[motherI] != 0) |
580 | { |
581 | rot = new G4RotationMatrix(*rot); |
582 | rot->rotateZ(-fCurrentPhiCenter[motherI]); |
583 | origin.rotateZ(fCurrentPhiCenter[motherI]); |
584 | } |
585 | |
586 | std::map<vpair_t,G4LogicalVolume*>::iterator mine; |
587 | for (mine = fLogicalVolumes.find(vpair_t(myvoluI,0)); |
588 | mine != fLogicalVolumes.end() && mine->first.first == myvoluI; |
589 | ++mine) |
590 | { |
591 | int ilayer = mine->first.second; |
592 | std::map<vpair_t,G4LogicalVolume*>::iterator moms; |
593 | moms = addNewLayer(motherI,fRef.fRelativeLayer + ilayer); |
594 | G4PVPlacement *pvol = new G4PVPlacement(rot, origin, |
595 | mine->second, |
596 | mine->second->GetName(), |
597 | moms->second,0,icopy); |
598 | #ifdef CHECK_OVERLAPS_MM |
599 | pvol->CheckOverlaps(1000,CHECK_OVERLAPS_MM); |
600 | #endif |
601 | #ifdef DEBUG_PLACEMENT |
602 | XString nameS(el->getAttribute(X("name")XString("name").unicode_str())); |
603 | G4cout(*G4cout_p) << "volume " << nameS |
604 | << "->" << mine->second->GetName() |
605 | << "->" << mine->second->GetSolid()->GetName() |
606 | << " being placed in mother " << moms->second->GetName() |
607 | << "->" << moms->second->GetSolid()->GetName() << " at " |
608 | << fRef.fOrigin[0]*cm << "," |
609 | << fRef.fOrigin[1]*cm << "," |
610 | << fRef.fOrigin[2]*cm |
611 | << G4endlstd::endl; |
612 | #endif |
613 | |
614 | if (ilayer == 0) { |
615 | vpair_t mycopy(myvoluI,icopy); |
616 | fPhysicalVolumes[mycopy] = pvol; |
617 | fCurrentPlacement[myvoluI] = fPhysicalVolumes.find(mycopy); |
618 | fCurrentMother[myvoluI] = moms; |
619 | } |
620 | } |
621 | |
622 | fPending = false; |
623 | } |
624 | |
625 | #ifdef LINUX_CPUTIME_PROFILING |
626 | timestr << " ( " << timer.getUserDelta() << " ) "; |
627 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
628 | #endif |
629 | return icopy; |
630 | } |
631 | |
632 | std::map<HddsG4Builder::vpair_t,G4LogicalVolume*>::iterator |
633 | HddsG4Builder::addNewLayer(int volume_id, int layer) |
634 | { |
635 | #ifdef LINUX_CPUTIME_PROFILING |
636 | std::stringstream timestr; |
637 | timestr << "addNewLayer: " << timer.getUserTime() << ", " |
638 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
639 | timer.resetClocks(); |
640 | #endif |
641 | |
642 | |
643 | |
644 | |
645 | |
646 | int max_layer = -1; |
647 | std::map<vpair_t,G4LogicalVolume*>::iterator volex; |
648 | for (volex = fLogicalVolumes.find(vpair_t(volume_id,0)); |
649 | volex != fLogicalVolumes.end() && |
650 | volex->first.first == volume_id; |
651 | ++volex) |
652 | { |
653 | max_layer = volex->first.second; |
654 | if (max_layer == layer) |
655 | { |
656 | return volex; |
657 | } |
658 | } |
659 | if (max_layer == -1) |
660 | { |
661 | G4cerr(*G4cerr_p) << "HddsG4Builder::addNewLayer called for volume " |
662 | << volume_id << " which does not (yet) exist!" |
663 | << G4endlstd::endl |
664 | << "Unable to continue, giving up." |
665 | << G4endlstd::endl; |
666 | exit(1); |
667 | } |
668 | |
669 | |
670 | |
671 | |
672 | |
673 | |
674 | |
675 | std::map<vpair_t,G4LogicalVolume*>::iterator moms; |
676 | if (fCurrentMother.find(volume_id) != fCurrentMother.end()) |
677 | { |
678 | moms = fCurrentMother[volume_id]; |
679 | moms = addNewLayer(moms->first.first, moms->first.second + layer); |
680 | } |
681 | else |
682 | { |
683 | moms = fLogicalVolumes.end(); |
684 | } |
685 | |
686 | |
687 | |
688 | |
689 | |
690 | |
691 | vpair_t newvol(volume_id,layer); |
692 | vpair_t oldvol(volume_id,max_layer); |
693 | G4VSolid *solid = fLogicalVolumes[oldvol]->GetSolid()->Clone(); |
694 | G4String name(fLogicalVolumes[vpair_t(volume_id,0)]->GetName()); |
695 | std::stringstream str; |
696 | str << name << "::" << layer; |
697 | solid->SetName(str.str()); |
698 | G4Material *material = fLogicalVolumes[oldvol]->GetMaterial(); |
699 | G4FieldManager *fieldmgr = fLogicalVolumes[oldvol]->GetFieldManager(); |
700 | if (layer > max_layer) |
701 | { |
702 | fLogicalVolumes[newvol] = new G4LogicalVolume(solid,material, |
703 | str.str(),fieldmgr); |
704 | fLogicalVolumes[oldvol]->SetMaterial(0); |
705 | } |
706 | else |
707 | { |
708 | fLogicalVolumes[newvol] = new G4LogicalVolume(solid,0,str.str(),fieldmgr); |
709 | } |
710 | |
711 | |
712 | |
713 | |
714 | if (moms == fLogicalVolumes.end()) |
715 | { |
716 | return fLogicalVolumes.find(newvol); |
717 | } |
718 | |
719 | |
720 | |
721 | |
722 | |
723 | |
724 | |
725 | |
726 | if (fCurrentPlacement.find(volume_id) != fCurrentPlacement.end()) |
727 | { |
728 | std::map<vpair_t,G4VPhysicalVolume*>::iterator placement; |
729 | placement = fCurrentPlacement[volume_id]; |
730 | G4PVPlacement *player0 = (G4PVPlacement*)placement->second; |
731 | vpair_t mycopy(volume_id,placement->first.second); |
732 | G4PVPlacement *playerN; |
733 | playerN = new G4PVPlacement(player0->GetRotation(), |
734 | player0->GetTranslation(), |
735 | fLogicalVolumes[newvol], |
736 | str.str(), |
737 | moms->second,0, |
738 | placement->first.second); |
739 | #ifdef CHECK_OVERLAPS_MM |
740 | playerN->CheckOverlaps(1000,CHECK_OVERLAPS_MM); |
741 | #endif |
742 | #ifdef DEBUG_PLACEMENT |
743 | G4cout(*G4cout_p) << "volume " << str.str() |
744 | << "->" << fLogicalVolumes[newvol]->GetName() |
745 | << "->" << fLogicalVolumes[newvol]->GetSolid()->GetName() |
746 | << " being placed in mother " << moms->second->GetName() |
747 | << "->" << moms->second->GetSolid()->GetName() |
748 | << G4endlstd::endl; |
749 | #endif |
750 | } |
751 | |
752 | |
753 | |
754 | if (fCurrentDivision.find(volume_id) != fCurrentDivision.end()) |
755 | { |
756 | std::map<vpair_t,G4VPhysicalVolume*>::iterator division; |
757 | division = fCurrentDivision[volume_id]; |
758 | G4PVDivision *division0 = (G4PVDivision*)division->second; |
759 | EAxis axis; |
760 | G4int ndiv; |
761 | G4double width; |
762 | G4double offset; |
763 | G4bool consuming; |
764 | division0->GetReplicationData(axis,ndiv,width,offset,consuming); |
765 | |
766 | axis = division0->GetDivisionAxis(); |
767 | vpair_t mydiv(volume_id,division->first.second); |
768 | G4PVDivision *divisionN; |
769 | divisionN = new G4PVDivision(str.str(), |
| Value stored to 'divisionN' is never read |
770 | fLogicalVolumes[newvol], |
771 | moms->second, |
772 | axis, ndiv, width, offset); |
773 | fLogicalVolumes[newvol]->SetVisAttributes(new G4VisAttributes(false)); |
774 | #ifdef DEBUG_PLACEMENT |
775 | G4cout(*G4cout_p) << ndiv << " copies of division " << str.str() |
776 | << "->" << fLogicalVolumes[newvol]->GetName() |
777 | << "->" << fLogicalVolumes[newvol]->GetSolid()->GetName() |
778 | << " being placed in mother " << moms->second->GetName() |
779 | << "->" << moms->second->GetSolid()->GetName() |
780 | << " offset by " << offset << " on axis " << axis |
781 | << " and repeating every " << width |
782 | << G4endlstd::endl; |
783 | #endif |
784 | } |
785 | |
786 | #ifdef LINUX_CPUTIME_PROFILING |
787 | timestr << " ( " << timer.getUserDelta() << " ) "; |
788 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
789 | #endif |
790 | return fLogicalVolumes.find(newvol); |
791 | } |
792 | |
793 | int HddsG4Builder::createDivision(XString& divStr, Refsys& ref) |
794 | { |
795 | #ifdef LINUX_CPUTIME_PROFILING |
796 | std::stringstream timestr; |
797 | timestr << "createDivision: " << timer.getUserTime() << ", " |
798 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
799 | timer.resetClocks(); |
800 | #endif |
801 | int ndiv = CodeWriter::createDivision(divStr,ref); |
802 | |
803 | XString motherS(ref.fMother->getAttribute(X("HDDSvolu")XString("HDDSvolu").unicode_str())); |
804 | XString myvoluS(ref.fPartition.divEl->getAttribute(X("HDDSvolu")XString("HDDSvolu").unicode_str())); |
805 | int motherI = atoi(S(motherS)motherS.c_str()); |
806 | int myvoluI = atoi(S(myvoluS)myvoluS.c_str()); |
807 | std::map<vpair_t,G4LogicalVolume*>::iterator moms; |
808 | moms = fLogicalVolumes.find(vpair_t(motherI,0)); |
809 | G4VSolid* solid = moms->second->GetSolid()->Clone(); |
810 | solid->SetName(divStr); |
811 | EAxis axis; |
812 | double width,offset; |
813 | |
814 | if (ref.fPartition.axis == "x") |
815 | { |
816 | axis = kXAxis; |
817 | width = ref.fPartition.step * cm; |
818 | offset = ref.fPartition.offset * cm; |
819 | } |
820 | else if (ref.fPartition.axis == "y") |
821 | { |
822 | axis = kYAxis; |
823 | width = ref.fPartition.step * cm; |
824 | offset = ref.fPartition.offset * cm; |
825 | } |
826 | else if (ref.fPartition.axis == "z") |
827 | { |
828 | axis = kZAxis; |
829 | width = ref.fPartition.step * cm; |
830 | offset = ref.fPartition.offset * cm; |
831 | } |
832 | else if (ref.fPartition.axis == "rho") |
833 | { |
834 | axis = kRho; |
835 | width = ref.fPartition.step * cm; |
836 | offset = ref.fPartition.offset * cm; |
837 | } |
838 | else if (ref.fPartition.axis == "phi") |
839 | { |
840 | axis = kPhi; |
841 | width = ref.fPartition.step * deg; |
842 | offset = ref.fPartition.offset * deg; |
843 | } |
844 | else |
845 | { |
846 | XString nameS(ref.fMother->getAttribute(X("name")XString("name").unicode_str())); |
847 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: Volume division is requested " |
848 | << " for volume " << S(nameS)nameS.c_str() << " but division axis " |
849 | << "\"" << ref.fPartition.axis << "\"" |
850 | << " is not currently supported." |
851 | << G4endlstd::endl; |
852 | exit(1); |
853 | } |
854 | |
855 | vpair_t mothvol(moms->first); |
856 | vpair_t myvol(myvoluI,moms->first.second); |
857 | G4Material* material = moms->second->GetMaterial(); |
858 | G4FieldManager* fieldmgr = fFieldManagers[ref.fRegionID]; |
859 | fLogicalVolumes[myvol] = new G4LogicalVolume(solid,material, |
860 | S(divStr)divStr.c_str(),fieldmgr); |
861 | fLogicalVolumes[myvol]->SetVisAttributes(new G4VisAttributes(false)); |
862 | vpair_t mydiv(myvoluI,0); |
863 | fPhysicalVolumes[mydiv] = new G4PVDivision(S(divStr)divStr.c_str(), |
864 | fLogicalVolumes[myvol], |
865 | fLogicalVolumes[mothvol], |
866 | axis,ndiv,width,offset); |
867 | fCurrentDivision[myvoluI] = fPhysicalVolumes.find(mydiv); |
868 | fCurrentMother[myvoluI] = moms; |
869 | |
870 | |
871 | |
872 | G4VPVParameterisation *param = fPhysicalVolumes[mydiv] |
873 | ->GetParameterisation(); |
874 | if (dynamic_cast<G4Box*>(solid)) |
875 | { |
876 | param->ComputeDimensions(*(G4Box*)solid,0,fPhysicalVolumes[mydiv]); |
877 | } |
878 | else if (dynamic_cast<G4Tubs*>(solid)) |
879 | { |
880 | param->ComputeDimensions(*(G4Tubs*)solid,0,fPhysicalVolumes[mydiv]); |
881 | double start = ((G4Tubs*)solid)->GetStartPhiAngle(); |
882 | double delta = ((G4Tubs*)solid)->GetDeltaPhiAngle(); |
883 | fCurrentPhiCenter[myvoluI] = start + delta/2; |
884 | } |
885 | else if (dynamic_cast<G4Trd*>(solid)) |
886 | { |
887 | param->ComputeDimensions(*(G4Trd*)solid,0,fPhysicalVolumes[mydiv]); |
888 | } |
889 | else if (dynamic_cast<G4Trap*>(solid)) |
890 | { |
891 | param->ComputeDimensions(*(G4Trap*)solid,0,fPhysicalVolumes[mydiv]); |
892 | } |
893 | else if (dynamic_cast<G4Cons*>(solid)) |
894 | { |
895 | param->ComputeDimensions(*(G4Cons*)solid,0,fPhysicalVolumes[mydiv]); |
896 | double start = ((G4Cons*)solid)->GetStartPhiAngle(); |
897 | double delta = ((G4Cons*)solid)->GetDeltaPhiAngle(); |
898 | fCurrentPhiCenter[myvoluI] = start + delta/2; |
899 | } |
900 | else if (dynamic_cast<G4Sphere*>(solid)) |
901 | { |
902 | param->ComputeDimensions(*(G4Sphere*)solid,0,fPhysicalVolumes[mydiv]); |
903 | double start = ((G4Sphere*)solid)->GetStartPhiAngle(); |
904 | double delta = ((G4Sphere*)solid)->GetDeltaPhiAngle(); |
905 | fCurrentPhiCenter[myvoluI] = start + delta/2; |
906 | } |
907 | else if (dynamic_cast<G4Orb*>(solid)) |
908 | { |
909 | param->ComputeDimensions(*(G4Orb*)solid,0,fPhysicalVolumes[mydiv]); |
910 | } |
911 | else if (dynamic_cast<G4Torus*>(solid)) |
912 | { |
913 | param->ComputeDimensions(*(G4Torus*)solid,0,fPhysicalVolumes[mydiv]); |
914 | double start = ((G4Torus*)solid)->GetSPhi(); |
915 | double delta = ((G4Torus*)solid)->GetDPhi(); |
916 | fCurrentPhiCenter[myvoluI] = start + delta/2; |
917 | } |
918 | else if (dynamic_cast<G4Para*>(solid)) |
919 | { |
920 | param->ComputeDimensions(*(G4Para*)solid,0,fPhysicalVolumes[mydiv]); |
921 | } |
922 | else if (dynamic_cast<G4Polycone*>(solid)) |
923 | { |
924 | param->ComputeDimensions(*(G4Polycone*)solid,0,fPhysicalVolumes[mydiv]); |
925 | G4PolyconeHistorical *params; |
926 | params = ((G4Polycone*)solid)->GetOriginalParameters(); |
927 | double start = params->Start_angle; |
928 | double delta = params->Opening_angle; |
929 | fCurrentPhiCenter[myvoluI] = start + delta/2; |
930 | } |
931 | else if (dynamic_cast<G4Polyhedra*>(solid)) |
932 | { |
933 | param->ComputeDimensions(*(G4Polyhedra*)solid,0,fPhysicalVolumes[mydiv]); |
934 | G4PolyhedraHistorical *params; |
935 | params = ((G4Polyhedra*)solid)->GetOriginalParameters(); |
936 | double start = params->Start_angle; |
937 | double delta = params->Opening_angle; |
938 | fCurrentPhiCenter[myvoluI] = start + delta/2; |
939 | } |
940 | else if (dynamic_cast<G4Hype*>(solid)) |
941 | { |
942 | param->ComputeDimensions(*(G4Hype*)solid,0,fPhysicalVolumes[mydiv]); |
943 | } |
944 | else |
945 | { |
946 | G4cerr(*G4cerr_p) << "HddsG4Builder::createDivision called for unsupported " |
947 | << "solid type!" << G4endlstd::endl |
948 | << "Unable to continue, giving up." |
949 | << G4endlstd::endl; |
950 | exit(1); |
951 | } |
952 | #ifdef DEBUG_PLACEMENT |
953 | G4cout(*G4cout_p) << ndiv << " copies of division " << divStr |
954 | << "->" << fLogicalVolumes[myvol]->GetName() |
955 | << "->" << fLogicalVolumes[myvol]->GetSolid()->GetName() |
956 | << " being placed in mother " << fLogicalVolumes[mothvol]->GetName() |
957 | << "->" << fLogicalVolumes[mothvol]->GetSolid()->GetName() |
958 | << " offset by " << offset << " on axis " << axis |
959 | << " and repeating every " << width |
960 | << G4endlstd::endl; |
961 | #endif |
962 | |
963 | #ifdef LINUX_CPUTIME_PROFILING |
964 | timestr << " ( " << timer.getUserDelta() << " ) "; |
965 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
966 | #endif |
967 | return ndiv; |
968 | } |
969 | |
970 | int HddsG4Builder::createRegion(DOMElement* el, Refsys& ref) |
971 | { |
972 | #ifdef LINUX_CPUTIME_PROFILING |
973 | std::stringstream timestr; |
974 | timestr << "createRegion: " << timer.getUserTime() << ", " |
975 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
976 | timer.resetClocks(); |
977 | #endif |
978 | int iregion = CodeWriter::createRegion(el,ref); |
979 | |
980 | XString motherS(ref.fMother->getAttribute(X("HDDSvolu")XString("HDDSvolu").unicode_str())); |
981 | int motherI = atoi(S(motherS)motherS.c_str()); |
982 | std::map<vpair_t,G4LogicalVolume*>::iterator moms; |
983 | moms = fLogicalVolumes.find(vpair_t(motherI,0)); |
984 | |
985 | G4double *R = (G4double*)ref.fMRmatrix; |
986 | G4double *O = ref.fMOrigin; |
987 | G4AffineTransform xform(CLHEP::HepRotation(CLHEP::HepRep3x3(R)), |
988 | CLHEP::Hep3Vector(O[0],O[1],O[2])); |
989 | |
990 | if (ref.fRegion) |
991 | { |
992 | DOMNodeList* noBfieldL; |
993 | DOMNodeList* uniBfieldL; |
994 | DOMNodeList* mapBfieldL; |
995 | DOMNodeList* compBfieldL; |
996 | DOMNodeList* swimL; |
997 | noBfieldL = ref.fRegion->getElementsByTagName(X("noBfield")XString("noBfield").unicode_str()); |
998 | uniBfieldL = ref.fRegion->getElementsByTagName(X("uniformBfield")XString("uniformBfield").unicode_str()); |
999 | mapBfieldL = ref.fRegion->getElementsByTagName(X("mappedBfield")XString("mappedBfield").unicode_str()); |
1000 | compBfieldL = ref.fRegion->getElementsByTagName(X("computedBfield")XString("computedBfield").unicode_str()); |
1001 | swimL = ref.fRegion->getElementsByTagName(X("swim")XString("swim").unicode_str()); |
1002 | double maxArcStep = 0; |
1003 | XString methodS("helix"); |
1004 | if (swimL->getLength() > 0) |
1005 | { |
1006 | DOMElement* swimEl = (DOMElement*)swimL->item(0); |
1007 | methodS = swimEl->getAttribute(X("method")XString("method").unicode_str()); |
1008 | Units unit; |
1009 | unit.getConversions(swimEl); |
1010 | XString stepS(swimEl->getAttribute(X("maxArcStep")XString("maxArcStep").unicode_str())); |
1011 | maxArcStep = atof(S(stepS)stepS.c_str())/unit.rad; |
1012 | } |
1013 | if (noBfieldL->getLength() > 0) |
1014 | { |
1015 | fMagneticRegions[iregion] = 0; |
1016 | G4ThreeVector Bvec(0,0,1e-99); |
1017 | G4UniformMagField *fld = new G4UniformMagField(Bvec); |
1018 | G4Mag_EqRhs *eqn = new G4Mag_UsualEqRhs(fld); |
1019 | G4MagIntegratorStepper *stepper = new G4ExactHelixStepper(eqn); |
1020 | G4ChordFinder *cfinder = new G4ChordFinder(fld, 1e+99, stepper); |
1021 | fFieldManagers[iregion] = new G4FieldManager(fld, cfinder); |
1022 | } |
1023 | else if (uniBfieldL->getLength() > 0) |
1024 | { |
1025 | Units funit; |
1026 | DOMElement* uniBfieldEl = (DOMElement*)uniBfieldL->item(0); |
1027 | funit.getConversions(uniBfieldEl); |
1028 | XString bvecS(uniBfieldEl->getAttribute(X("Bx_By_Bz")XString("Bx_By_Bz").unicode_str())); |
1029 | std::stringstream str(S(bvecS)bvecS.c_str()); |
1030 | double B[3]; |
1031 | str >> B[0] >> B[1] >> B[2]; |
1032 | double u = kilogauss/funit.kG; |
1033 | G4ThreeVector Bvec(B[0],B[1],B[2]); |
1034 | GlueXUniformMagField *fld = new GlueXUniformMagField(Bvec,u,xform); |
1035 | fMagneticRegions[iregion] = fld; |
1036 | G4Mag_EqRhs *eqn = new G4Mag_UsualEqRhs(fld); |
1037 | G4MagIntegratorStepper *stepper = new G4ExactHelixStepper(eqn); |
1038 | G4ChordFinder *cfinder = new G4ChordFinder(fld, 0.01, stepper); |
1039 | fFieldManagers[iregion] = new G4FieldManager(fld, cfinder); |
1040 | } |
1041 | else if (mapBfieldL->getLength() > 0) |
1042 | { |
1043 | Units funit; |
1044 | DOMElement* mapBfieldEl = (DOMElement*)mapBfieldL->item(0); |
1045 | funit.getConversions(mapBfieldEl); |
1046 | XString bvecS(mapBfieldEl->getAttribute(X("maxBfield")XString("maxBfield").unicode_str())); |
1047 | std::stringstream str(S(bvecS)bvecS.c_str()); |
1048 | double Bmax; |
1049 | str >> Bmax; |
1050 | double u = kilogauss/funit.kG; |
1051 | GlueXMappedMagField *fld = new GlueXMappedMagField(Bmax,u,xform); |
1052 | fMagneticRegions[iregion] = fld; |
1053 | G4Mag_EqRhs *eqn = new G4Mag_UsualEqRhs(fld); |
1054 | G4MagIntegratorStepper *stepper; |
1055 | if (methodS == "RungeKutta") { |
1056 | stepper = new G4ClassicalRK4(eqn); |
1057 | } |
1058 | else { |
1059 | stepper = new G4HelixMixedStepper(eqn); |
1060 | } |
1061 | G4ChordFinder *cfinder = new G4ChordFinder(fld, 0.01, stepper); |
1062 | if (maxArcStep > 0) { |
1063 | double rmin = 0.1 / (0.03 * Bmax * u) * meter; |
1064 | double max_miss = rmin * (1 - cos(maxArcStep / 2)); |
1065 | cfinder->SetDeltaChord(max_miss); |
1066 | } |
1067 | fFieldManagers[iregion] = new G4FieldManager(fld, cfinder); |
1068 | } |
1069 | else if (compBfieldL->getLength() > 0) |
1070 | { |
1071 | Units funit; |
1072 | DOMElement* compBfieldEl = (DOMElement*)compBfieldL->item(0); |
1073 | funit.getConversions(compBfieldEl); |
1074 | XString bvecS(compBfieldEl->getAttribute(X("maxBfield")XString("maxBfield").unicode_str())); |
1075 | std::stringstream str(S(bvecS)bvecS.c_str()); |
1076 | double Bmax; |
1077 | str >> Bmax; |
1078 | double u = kilogauss/funit.kG; |
1079 | GlueXComputedMagField *fld = new GlueXComputedMagField(Bmax,u,xform); |
1080 | fld->SetFunction(XString(compBfieldEl->getAttribute(X("function")XString("function").unicode_str()))); |
1081 | fMagneticRegions[iregion] = fld; |
1082 | G4Mag_EqRhs *eqn = new G4Mag_UsualEqRhs(fld); |
1083 | G4MagIntegratorStepper *stepper; |
1084 | if (methodS == "RungeKutta") { |
1085 | stepper = new G4ClassicalRK4(eqn); |
1086 | } |
1087 | else { |
1088 | stepper = new G4HelixMixedStepper(eqn); |
1089 | } |
1090 | G4ChordFinder *cfinder = new G4ChordFinder(fld, 0.01, stepper); |
1091 | if (maxArcStep > 0) { |
1092 | double rmin = 0.1 / (0.03 * Bmax * u) * meter; |
1093 | double max_miss = rmin * (1 - cos(maxArcStep / 2)); |
1094 | cfinder->SetDeltaChord(max_miss); |
1095 | } |
1096 | fFieldManagers[iregion] = new G4FieldManager(fld, cfinder); |
1097 | } |
1098 | } |
1099 | |
1100 | #ifdef LINUX_CPUTIME_PROFILING |
1101 | timestr << " ( " << timer.getUserDelta() << " ) "; |
1102 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
1103 | #endif |
1104 | return iregion; |
1105 | } |
1106 | |
1107 | void HddsG4Builder::createSetFunctions(DOMElement* el, const XString& ident) |
1108 | { |
1109 | #ifdef LINUX_CPUTIME_PROFILING |
1110 | std::stringstream timestr; |
1111 | timestr << "createSetFunctions: " << timer.getUserTime() << ", " |
1112 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
1113 | timer.resetClocks(); |
1114 | #endif |
1115 | CodeWriter::createSetFunctions(el,ident); |
1116 | |
1117 | #ifdef LINUX_CPUTIME_PROFILING |
1118 | timestr << " ( " << timer.getUserDelta() << " ) "; |
1119 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
1120 | #endif |
1121 | } |
1122 | |
1123 | void HddsG4Builder::createGetFunctions(DOMElement* el, const XString& ident) |
1124 | { |
1125 | #ifdef LINUX_CPUTIME_PROFILING |
1126 | std::stringstream timestr; |
1127 | timestr << "createGetFunctions: " << timer.getUserTime() << ", " |
1128 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
1129 | timer.resetClocks(); |
1130 | #endif |
1131 | CodeWriter::createGetFunctions(el,ident); |
1132 | |
1133 | std::vector<int> table; |
1134 | std::vector<int> start; |
1135 | start.push_back(0); |
1136 | |
1137 | XString funcNameStr; |
1138 | XString identCaps(ident); |
1139 | identCaps[0] = toupper(identCaps[0]); |
1140 | funcNameStr = "get" + identCaps; |
1141 | for (int ivolu = 1; ivolu <= Refsys::fVolumes; ivolu++) |
1142 | { |
1143 | start.push_back(0); |
1144 | int ncopy = Refsys::fIdentifierTable[ivolu]["copy counter"].back(); |
1145 | std::map<std::string,std::vector<int> >::iterator idlist = |
1146 | Refsys::fIdentifierTable[ivolu].find(ident); |
1147 | if (idlist != Refsys::fIdentifierTable[ivolu].end()) |
1148 | { |
1149 | if (ncopy != (int)idlist->second.size()) |
1150 | { |
1151 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " warning: volume " << ivolu |
1152 | << " has " << ncopy << " copies, but " |
1153 | << idlist->second.size() << " " |
1154 | << ident << " identifiers!" << G4endlstd::endl; |
1155 | for (int idx = 0; idx < (int)idlist->second.size(); idx++) |
1156 | { |
1157 | G4cerr(*G4cerr_p) << idlist->second[idx] << " "; |
1158 | if (idx/20*20 == idx) |
1159 | G4cerr(*G4cerr_p) << G4endlstd::endl; |
1160 | } |
1161 | G4cerr(*G4cerr_p) << G4endlstd::endl; |
1162 | } |
1163 | start[ivolu] = table.size() + 1; |
1164 | for (int icopy = 0; icopy < ncopy; icopy++) |
1165 | { |
1166 | table.push_back(idlist->second[icopy]); |
1167 | } |
1168 | } |
1169 | else |
1170 | { |
1171 | start[ivolu] = 0; |
1172 | } |
1173 | } |
1174 | |
1175 | #ifdef LINUX_CPUTIME_PROFILING |
1176 | timestr << " ( " << timer.getUserDelta() << " ) "; |
1177 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
1178 | #endif |
1179 | } |
1180 | |
1181 | void HddsG4Builder::createMapFunctions(DOMElement* el, const XString& ident) |
1182 | { |
1183 | #ifdef LINUX_CPUTIME_PROFILING |
1184 | std::stringstream timestr; |
1185 | timestr << "createMapFunctions: " << timer.getUserTime() << ", " |
1186 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
1187 | timer.resetClocks(); |
1188 | #endif |
1189 | CodeWriter::createMapFunctions(el,ident); |
1190 | |
1191 | if (el == 0) |
1192 | { |
1193 | return; |
1194 | } |
1195 | DOMNodeList* regionsL = el->getOwnerDocument()->getDocumentElement() |
1196 | ->getElementsByTagName(X("regions")XString("regions").unicode_str()); |
1197 | if (regionsL->getLength() == 0) |
1198 | { |
1199 | return; |
1200 | } |
1201 | DOMElement* regionsEl = (DOMElement*)regionsL->item(0); |
1202 | DOMNodeList* regionL = regionsEl->getElementsByTagName(X("region")XString("region").unicode_str()); |
1203 | for (int ireg=0; ireg < (int)regionL->getLength(); ++ireg) |
1204 | { |
1205 | DOMElement* regionEl = (DOMElement*)regionL->item(ireg); |
1206 | DOMNodeList* mapfTagL = regionEl->getElementsByTagName(X("mappedBfield")XString("mappedBfield").unicode_str()); |
1207 | if (mapfTagL->getLength() == 0) |
1208 | continue; |
1209 | DOMElement* mapfEl = (DOMElement*)mapfTagL->item(0); |
1210 | XString nameS(regionEl->getAttribute(X("name")XString("name").unicode_str())); |
1211 | XString iregionS(regionEl->getAttribute(X("HDDSregion")XString("HDDSregion").unicode_str())); |
1212 | if (iregionS.size() == 0) |
1213 | continue; |
1214 | int iregion = atoi(S(iregionS)iregionS.c_str()); |
1215 | GlueXMappedMagField *magfield = (GlueXMappedMagField*) |
1216 | fMagneticRegions[iregion]; |
1217 | int axorder[] = {0,0,0,0}; |
1218 | int axsamples[] = {0,0,0,0}; |
1219 | |
1220 | XString gridtype; |
1221 | DOMNodeList* gridL = mapfEl->getElementsByTagName(X("grid")XString("grid").unicode_str()); |
1222 | int ngrid; |
1223 | for (ngrid = 0; ngrid < (int)gridL->getLength(); ++ngrid) |
1224 | { |
1225 | int axsense[] = {1,1,1,1}; |
1226 | double axlower[4], axupper[4]; |
1227 | DOMElement* gridEl = (DOMElement*)gridL->item(ngrid); |
1228 | XString typeS(gridEl->getAttribute(X("type")XString("type").unicode_str())); |
1229 | if (gridtype.size() > 0 && typeS != gridtype) |
1230 | { |
1231 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: mappedBfield in region " << S(nameS)nameS.c_str() |
1232 | << " superimposes incompatible grid types." << G4endlstd::endl; |
1233 | exit(1); |
1234 | } |
1235 | gridtype = typeS; |
1236 | |
1237 | DOMNodeList* samplesL = gridEl->getElementsByTagName(X("samples")XString("samples").unicode_str()); |
1238 | if (samplesL->getLength() != 3) |
1239 | { |
1240 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: mappedBfield in region " << S(nameS)nameS.c_str() |
1241 | << " does not have samples for three axes." << G4endlstd::endl; |
1242 | exit(1); |
1243 | } |
1244 | |
1245 | for (int iax = 1; iax <= 3; ++iax) |
1246 | { |
1247 | DOMElement* sampleEl = (DOMElement*)samplesL->item(iax-1); |
1248 | XString nS(sampleEl->getAttribute(X("n")XString("n").unicode_str())); |
1249 | XString axisS(sampleEl->getAttribute(X("axis")XString("axis").unicode_str())); |
1250 | XString boundsS(sampleEl->getAttribute(X("bounds")XString("bounds").unicode_str())); |
1251 | XString senseS(sampleEl->getAttribute(X("sense")XString("sense").unicode_str())); |
1252 | Units sunit; |
1253 | double bound[2]; |
1254 | sunit.getConversions(sampleEl); |
1255 | std::stringstream listr(boundsS); |
1256 | listr >> bound[0] >> bound[1]; |
1257 | int iaxis=0; |
1258 | if (gridtype == "cartesian") |
1259 | { |
1260 | if (axisS == "x" && |
1261 | (axorder[0] == 0 || axorder[0] == iax)) |
1262 | { |
1263 | iaxis = 1; |
1264 | axorder[0] = iax; |
1265 | bound[0] /= sunit.cm; |
1266 | bound[1] /= sunit.cm; |
1267 | } |
1268 | else if (axisS == "y" && |
1269 | (axorder[1] == 0 || axorder[1] == iax)) |
1270 | { |
1271 | iaxis = 2; |
1272 | axorder[1] = iax; |
1273 | bound[0] /= sunit.cm; |
1274 | bound[1] /= sunit.cm; |
1275 | } |
1276 | else if (axisS == "z" && |
1277 | (axorder[2] == 0 || axorder[2] == iax)) |
1278 | { |
1279 | iaxis = 3; |
1280 | axorder[2] = iax; |
1281 | bound[0] /= sunit.cm; |
1282 | bound[1] /= sunit.cm; |
1283 | } |
1284 | else |
1285 | { |
1286 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: grid in region " << S(nameS)nameS.c_str() |
1287 | << " contains an incompatible set of samples." |
1288 | << G4endlstd::endl; |
1289 | exit(1); |
1290 | } |
1291 | } |
1292 | else if (gridtype == "cylindrical") |
1293 | { |
1294 | if (axisS == "r" && axorder[0] == 0) |
1295 | { |
1296 | iaxis = 1; |
1297 | axorder[0] = iax; |
1298 | bound[0] /= sunit.cm; |
1299 | bound[1] /= sunit.cm; |
1300 | } |
1301 | else if (axisS == "phi" && axorder[1] == 0) |
1302 | { |
1303 | iaxis = 2; |
1304 | axorder[1] = iax; |
1305 | bound[0] /= sunit.rad; |
1306 | bound[1] /= sunit.rad; |
1307 | } |
1308 | else if (axisS == "z" && axorder[2] == 0) |
1309 | { |
1310 | iaxis = 3; |
1311 | axorder[2] = iax; |
1312 | bound[0] /= sunit.cm; |
1313 | bound[1] /= sunit.cm; |
1314 | } |
1315 | else |
1316 | { |
1317 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: grid in region " << S(nameS)nameS.c_str() |
1318 | << " contains an incompatible set of samples." |
1319 | << G4endlstd::endl; |
1320 | exit(1); |
1321 | } |
1322 | } |
1323 | |
1324 | int n = atoi(S(nS)nS.c_str()); |
1325 | if (axsamples[iaxis] == 0 || axsamples[iaxis] == n) |
1326 | { |
1327 | axsamples[iaxis] = n; |
1328 | } |
1329 | else |
1330 | { |
1331 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: mappedBfield in region " << S(nameS)nameS.c_str() |
1332 | << " combines incompatible grid elements." |
1333 | << G4endlstd::endl; |
1334 | exit(1); |
1335 | } |
1336 | |
1337 | axlower[iaxis] = bound[0]; |
1338 | axupper[iaxis] = bound[1]; |
1339 | if (senseS == "reverse") |
1340 | { |
1341 | axsense[iaxis] = -1; |
1342 | } |
1343 | } |
1344 | |
1345 | if (gridtype == "cartesian") |
1346 | { |
1347 | magfield->AddCartesianGrid(axsamples, axorder, axsense, |
1348 | axlower, axupper); |
1349 | } |
1350 | else if (gridtype == "cylindrical") |
1351 | { |
1352 | magfield->AddCylindricalGrid(axsamples, axorder, axsense, |
1353 | axlower, axupper); |
1354 | } |
1355 | else |
1356 | { |
1357 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: unrecognized grid type " |
1358 | << S(gridtype)gridtype.c_str() << G4endlstd::endl; |
1359 | exit(1); |
1360 | } |
1361 | } |
1362 | |
1363 | XString mapS(mapfEl->getAttribute(X("map")XString("map").unicode_str())); |
1364 | XString encS(mapfEl->getAttribute(X("encoding")XString("encoding").unicode_str())); |
1365 | if (encS != "utf-8") |
1366 | { |
1367 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: mappedBfield in region " << S(nameS)nameS.c_str() |
1368 | << " uses unsupported encoding " << encS << G4endlstd::endl; |
1369 | exit(1); |
1370 | } |
1371 | else if (mapS.substr(0,7) != "file://") |
1372 | { |
1373 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: mappedBfield in region " << S(nameS)nameS.c_str() |
1374 | << " uses unsupported map URL " << mapS << G4endlstd::endl; |
1375 | exit(1); |
1376 | } |
1377 | mapS.erase(0,7); |
1378 | magfield->ReadMapFile(S(mapS)mapS.c_str()); |
1379 | } |
1380 | |
1381 | |
1382 | |
1383 | |
1384 | addReflections(1); |
1385 | |
1386 | #ifdef LINUX_CPUTIME_PROFILING |
1387 | timestr << " ( " << timer.getUserDelta() << " ) "; |
1388 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
1389 | #endif |
1390 | } |
1391 | |
1392 | void HddsG4Builder::translate(DOMElement* topel) |
1393 | { |
1394 | #ifdef LINUX_CPUTIME_PROFILING |
1395 | std::stringstream timestr; |
1396 | timestr << "translate: " << timer.getUserTime() << ", " |
1397 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
1398 | timer.resetClocks(); |
1399 | #endif |
1400 | |
1401 | CodeWriter::translate(topel); |
1402 | DOMDocument* document = topel->getOwnerDocument(); |
1403 | DOMElement* worldEl = document->getElementById(X("everything")XString("everything").unicode_str()); |
1404 | XString topvolS(worldEl->getAttribute(X("envelope")XString("envelope").unicode_str())); |
1405 | DOMElement* topvolEl = document->getElementById(X(topvolS)XString(topvolS).unicode_str()); |
1406 | XString ivoluS(topvolEl->getAttribute(X("HDDSvolu")XString("HDDSvolu").unicode_str())); |
1407 | if (ivoluS.size() == 0) |
1408 | { |
1409 | G4cerr(*G4cerr_p) << APP_NAME"HddsG4Builder" << " error: top-level composition named " |
1410 | << "\"everything\" was not found in input document!" |
1411 | << G4endlstd::endl; |
1412 | exit(1); |
1413 | } |
1414 | fWorldVolume = atoi(S(ivoluS)ivoluS.c_str()); |
1415 | |
1416 | #ifdef LINUX_CPUTIME_PROFILING |
1417 | timestr << " ( " << timer.getUserDelta() << " ) "; |
1418 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
1419 | #endif |
1420 | } |
1421 | |
1422 | G4LogicalVolume* HddsG4Builder::getWorldVolume(int parallel) const |
1423 | { |
1424 | int worlds = 0; |
1425 | std::map<vpair_t,G4LogicalVolume*>::const_iterator paraworld; |
1426 | for (paraworld = fLogicalVolumes.find(vpair_t(fWorldVolume,0)); |
1427 | paraworld != fLogicalVolumes.end() && |
1428 | paraworld->first.first == fWorldVolume; |
1429 | ++paraworld, ++worlds) |
1430 | {} |
1431 | if (parallel < worlds) |
1432 | { |
1433 | for (int p=0; p <= parallel; ++p, --paraworld) {} |
1434 | return paraworld->second; |
1435 | } |
1436 | return 0; |
1437 | } |
1438 | |
1439 | void HddsG4Builder::addReflections(int volume_id) |
1440 | { |
1441 | #ifdef LINUX_CPUTIME_PROFILING |
1442 | std::stringstream timestr; |
1443 | timestr << "addReflections: " << timer.getUserTime() << ", " |
1444 | << timer.getSystemTime() << ", " << timer.getRealTime(); |
1445 | timer.resetClocks(); |
1446 | #endif |
1447 | std::vector<int> layer; |
1448 | std::map<vpair_t,G4LogicalVolume*>::iterator witer; |
1449 | for (witer = fLogicalVolumes.find(vpair_t(1,0)); |
1450 | witer != fLogicalVolumes.end() && witer->first.first == 1; |
1451 | ++witer) |
1452 | { |
1453 | layer.push_back(witer->first.second); |
1454 | } |
1455 | |
1456 | std::map<int,G4FieldManager*> lastManager; |
1457 | std::map<int,G4FieldManager*>::iterator miter; |
1458 | std::map<vpair_t,G4LogicalVolume*>::iterator mine; |
1459 | for (mine = fLogicalVolumes.find(vpair_t(volume_id,0)); |
1460 | mine != fLogicalVolumes.end() && mine->first.first == volume_id; |
1461 | ++mine) |
1462 | { |
1463 | G4FieldManager* fieldmgr = mine->second->GetFieldManager(); |
1464 | for (int child = 0; child < (int)mine->second->GetNoDaughters(); ++child) { |
1465 | G4VPhysicalVolume* pvol = mine->second->GetDaughter(child); |
1466 | G4LogicalVolume* lvol = pvol->GetLogicalVolume(); |
1467 | G4FieldManager* mgr = lvol->GetFieldManager(); |
1468 | int childId = getVolumeId(lvol); |
1469 | miter = lastManager.find(childId); |
1470 | if ((miter != lastManager.end() && mgr != miter->second) || |
1471 | (miter == lastManager.end() && mine->first.second != 0 && |
1472 | mgr != 0 && mgr != fieldmgr)) |
1473 | { |
1474 | G4cerr(*G4cerr_p) << "HddsG4Builder::addReflections warning - " |
1475 | << "local magnetic field assigned to volume " |
1476 | << lvol->GetName() |
1477 | << " placed on a buried layer in the geometry." |
1478 | << G4endlstd::endl |
1479 | << "This is a violation of the assumption made " |
1480 | << "in the Geant4 tracking algorithm that fields " |
1481 | << "at a given point are the same on all layers." |
1482 | << G4endlstd::endl |
1483 | << "Expect inconsistent results." |
1484 | << G4endlstd::endl; |
1485 | } |
1486 | else { |
1487 | lastManager[childId] = mgr; |
1488 | } |
1489 | |
1490 | |
1491 | |
1492 | if (miter == lastManager.end()) { |
1493 | addReflections(childId); |
1494 | } |
1495 | |
1496 | |
1497 | |
1498 | |
1499 | |
1500 | if (mgr != 0 && mgr != fieldmgr) { |
1501 | vpair_t cpy(childId, pvol->GetCopyNo()); |
1502 | std::map<vpair_t,G4VPhysicalVolume*>::iterator piter; |
1503 | piter = fPhysicalVolumes.find(cpy); |
1504 | if (piter == fPhysicalVolumes.end()) { |
1505 | G4cerr(*G4cerr_p) << "HddsG4Builder::addReflections error - " |
1506 | << "physical volume " << pvol->GetName() |
1507 | << " copy " << cpy.second |
1508 | << " not found in physical volume lookup table, " |
1509 | << "cannot continue." << G4endlstd::endl; |
1510 | exit(1); |
1511 | } |
1512 | else if (piter->second->IsReplicated()) { |
1513 | fCurrentDivision[childId] = piter; |
1514 | } |
1515 | else { |
1516 | fCurrentPlacement[childId] = piter; |
1517 | } |
1518 | fCurrentMother[childId] = mine; |
1519 | std::vector<int>::iterator liter; |
1520 | for (liter = layer.begin(); liter != layer.end(); ++liter) { |
1521 | if (*liter > mine->first.second) { |
1522 | #if DEBUG_REFLECTION |
1523 | G4cout(*G4cout_p) << "HddsG4Builder::addReflections info - " |
1524 | << "reflected volume " << lvol->GetName() |
1525 | << " copy " << cpy.second |
1526 | << " in mother " << mine->second->GetName() |
1527 | << " onto layer " << *liter |
1528 | << " because child field manager " << mgr |
1529 | << " is different from mother field manager " << fieldmgr |
1530 | << G4endlstd::endl; |
1531 | #endif |
1532 | addNewLayer(childId, *liter); |
1533 | } |
1534 | } |
1535 | } |
1536 | } |
1537 | } |
1538 | #ifdef LINUX_CPUTIME_PROFILING |
1539 | timestr << " ( " << timer.getUserDelta() << " ) "; |
1540 | G4cerr(*G4cerr_p) << timestr.str() << G4endlstd::endl; |
1541 | #endif |
1542 | } |
1543 | |
1544 | int HddsG4Builder::getVolumeId(G4LogicalVolume* vol) const |
1545 | { |
1546 | std::map<vpair_t,G4LogicalVolume*>::const_iterator iter; |
1547 | for (iter = fLogicalVolumes.begin(); iter != fLogicalVolumes.end(); ++iter) |
1548 | { |
1549 | if (iter->second == vol) { |
1550 | return iter->first.first; |
1551 | } |
1552 | } |
1553 | G4cerr(*G4cerr_p) << "HddsG4Builder::getVolumeId error - " |
1554 | << "logical volume " << vol->GetName() |
1555 | << " not found in volume lookup table, " |
1556 | << "cannot continue." << G4endlstd::endl; |
1557 | exit(1); |
1558 | } |
1559 | |
1560 | G4LogicalVolume* HddsG4Builder::getVolume(const G4String volname) const |
1561 | { |
1562 | std::map<vpair_t,G4LogicalVolume*>::const_iterator iter; |
1563 | for (iter = fLogicalVolumes.begin(); iter != fLogicalVolumes.end(); ++iter) |
1564 | { |
1565 | if (iter->second->GetName() == volname) { |
1566 | return iter->second; |
1567 | } |
1568 | } |
1569 | return 0; |
1570 | } |
1571 | |
1572 | const std::map<int, G4LogicalVolume*> HddsG4Builder::getSensitiveVolumes() const |
1573 | { |
1574 | return fSensitiveVolumes; |
1575 | } |