92 {
93
94 m_sensitive =
new CDCSensitiveDetector(
"CDCSensitiveDetector", (2 * 24)* CLHEP::eV, 10 * CLHEP::MeV);
95
96 const G4double realTemperture = (273.15 + 23.) * CLHEP::kelvin;
105
106 G4double h2odensity = 1.000 * CLHEP::g / CLHEP::cm3;
107 G4double a = 1.01 * CLHEP::g / CLHEP::mole;
108 G4Element* elH = new G4Element("Hydrogen", "H", 1., a);
109 a = 16.00 * CLHEP::g / CLHEP::mole;
110 G4Element* elO = new G4Element("Oxygen", "O", 8., a);
111 G4Material* medH2O = new G4Material("Water", h2odensity, 2);
112 medH2O->AddElement(elH, 2);
113 medH2O->AddElement(elO, 1);
116
117
118
119
120
121 const double rmax_innerWall = geo.getFiducialRmin();
122 const double rmin_outerWall = geo.getFiducialRmax();
123 const double diameter_senseWire = geo.getSenseDiameter();
124 const double diameter_fieldWire = geo.getFieldDiameter();
125 const double num_senseWire = static_cast<double>(geo.getNSenseWires());
126 const double num_fieldWire = static_cast<double>(geo.getNFieldWires());
127 double totalCS = M_PI * (rmin_outerWall * rmin_outerWall - rmax_innerWall * rmax_innerWall);
128
129
130 double senseCS = M_PI * (diameter_senseWire / 2) * (diameter_senseWire / 2) * num_senseWire;
131
132
133 double fieldCS = M_PI * (diameter_fieldWire / 2) * (diameter_fieldWire / 2) * num_fieldWire;
134
135
136 const double denHelium = medHelium->GetDensity() / 2.0;
137 const double denEthane = medEthane->GetDensity() / 2.0;
138 const double denAluminum = medAluminum->GetDensity() * (fieldCS / totalCS);
139 const double denTungsten = medTungsten->GetDensity() * (senseCS / totalCS);
140 const double density = denHelium + denEthane + denAluminum + denTungsten;
141 G4Material* cdcMed = new G4Material("CDCGasWire", density, 4, kStateGas, realTemperture);
142 cdcMed->AddMaterial(medHelium, denHelium / density);
143 cdcMed->AddMaterial(medEthane, denEthane / density);
144 cdcMed->AddMaterial(medTungsten, denTungsten / density);
145 cdcMed->AddMaterial(medAluminum, denAluminum / density);
146
147 G4Material* cdcMedGas = cdcMed;
148
151
152
153
154 if (gcp.getMaterialDefinitionMode() == 2) {
155 const double density2 = denHelium + denEthane;
156 cdcMedGas = new G4Material("CDCRealGas", density2, 2, kStateGas, realTemperture);
157 cdcMedGas->AddMaterial(medHelium, denHelium / density2);
158 cdcMedGas->AddMaterial(medEthane, denEthane / density2);
159 }
160
161 if (gcp.getPrintMaterialTable()) {
162 G4cout << *(G4Material::GetMaterialTable());
163 }
164
165 const auto& mother = geo.getMotherVolume();
166 const auto& motherRmin = mother.getRmin();
167 const auto& motherRmax = mother.getRmax();
168 const auto& motherZ = mother.getZ();
169 G4Polycone* solid_cdc =
170 new G4Polycone("solidCDC", 0 * CLHEP::deg, 360.* CLHEP::deg,
171 mother.getNNodes(), motherZ.data(),
172 motherRmin.data(), motherRmax.data());
173 m_logicalCDC =
new G4LogicalVolume(solid_cdc, medAir,
"logicalCDC", 0, 0, 0);
174 m_physicalCDC =
new G4PVPlacement(0, G4ThreeVector(geo.getGlobalOffsetX() * CLHEP::cm,
175 geo.getGlobalOffsetY() * CLHEP::cm,
177 "physicalCDC", &topVolume, false, 0);
178
179
180 G4Region* aRegion = new G4Region("CDCEnvelope");
183
184 m_VisAttributes.push_back(
new G4VisAttributes(
true, G4Colour(0., 1., 0.)));
185 for (const auto& wall : geo.getOuterWalls()) {
186 const int iOuterWall = wall.getId();
187 const string wallName = wall.getName();
188 const double wallRmin = wall.getRmin();
189 const double wallRmax = wall.getRmax();
190 const double wallZfwd = wall.getZfwd();
191 const double wallZbwd = wall.getZbwd();
192 const double length = (wallZfwd - wallZbwd) / 2.0;
193
194
195 G4Material* medWall;
196 if (strstr((wallName).c_str(), "MiddleWall") != nullptr) {
197 medWall = medCFRP;
198 } else {
199 medWall = medAluminum;
200 }
201 G4Tubs* outerWallTubeShape = new G4Tubs("solid" + wallName, wallRmin * CLHEP::cm,
202 wallRmax * CLHEP::cm, length * CLHEP::cm, 0 * CLHEP::deg, 360.*CLHEP::deg);
203
204 G4LogicalVolume* outerWallTube = new G4LogicalVolume(outerWallTubeShape, medWall, "solid" + wallName, 0, 0, 0);
206 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (length + wallZbwd)*CLHEP::cm), outerWallTube, "logical" + wallName,
208 }
209
210
211 m_VisAttributes.push_back(
new G4VisAttributes(
true, G4Colour(0., 1., 0.)));
212 for (const auto& wall : geo.getInnerWalls()) {
213 const string wallName = wall.getName();
214 const double wallRmin = wall.getRmin();
215 const double wallRmax = wall.getRmax();
216 const double wallZfwd = wall.getZfwd();
217 const double wallZbwd = wall.getZbwd();
218 const double length = (wallZfwd - wallZbwd) / 2.0;
219 const int iInnerWall = wall.getId();
220
221 G4Material* medWall;
222 if (strstr(wallName.c_str(), "MiddleWall") != nullptr) {
223 medWall = medCFRP;
224 } else if (strstr(wallName.c_str(), "MiddleGlue") != nullptr) {
225 medWall = medGlue;
226 } else {
227 medWall = medAluminum;
228 }
229
230 G4Tubs* innerWallTubeShape = new G4Tubs("solid" + wallName, wallRmin * CLHEP::cm,
231 wallRmax * CLHEP::cm, length * CLHEP::cm, 0 * CLHEP::deg, 360.*CLHEP::deg);
232 G4LogicalVolume* innerWallTube = new G4LogicalVolume(innerWallTubeShape, medWall, "logical" + wallName, 0, 0, 0);
234 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (length + wallZbwd)*CLHEP::cm), innerWallTube, "physical" + wallName,
236
237
238 }
239
240
241
242
243
244
245 const uint nSLayer = geo.getNSenseLayers();
246 const double length_feedthrough = geo.getFeedthroughLength();
247 for (uint iSLayer = 0; iSLayer < nSLayer; ++iSLayer) {
248 const auto& endplate = geo.getEndPlate(iSLayer);
249 const int nEPLayer = endplate.getNEndPlateLayers();
250
251 double rmin_sensitive_left, rmax_sensitive_left;
252 double rmin_sensitive_middle, rmax_sensitive_middle;
253 double rmin_sensitive_right, rmax_sensitive_right;
254 double zback_sensitive_left, zfor_sensitive_left;
255 double zback_sensitive_middle, zfor_sensitive_middle;
256 double zback_sensitive_right, zfor_sensitive_right;
257
258 if (iSLayer == 0) {
259 const auto& epLayerBwd = endplate.getEndPlateLayer(1);
260 const auto& epLayerFwd = endplate.getEndPlateLayer((nEPLayer / 2) + 1);
261 const auto& senseLayer = geo.getSenseLayer(iSLayer);
262 const auto& fieldLayer = geo.getFieldLayer(iSLayer);
263
264 rmin_sensitive_left = epLayerBwd.getRmax();
265 rmax_sensitive_left = fieldLayer.getR();
266 zback_sensitive_left = senseLayer.getZbwd();
267 zfor_sensitive_left = epLayerBwd.getZfwd();
268
269 rmin_sensitive_middle = (geo.getInnerWall(0)).getRmax();
270 rmax_sensitive_middle = fieldLayer.getR();
271 zback_sensitive_middle = epLayerBwd.getZfwd();
272 zfor_sensitive_middle = epLayerFwd.getZbwd();
273
274 rmin_sensitive_right = epLayerFwd.getRmax();
275 rmax_sensitive_right = fieldLayer.getR();
276 zback_sensitive_right = epLayerFwd.getZbwd();
277 zfor_sensitive_right = senseLayer.getZfwd();
278 } else if (iSLayer >= 1 && iSLayer <= 6) {
279 const auto& epLayerBwd = endplate.getEndPlateLayer(1);
280 const auto& epLayerFwd = endplate.getEndPlateLayer((nEPLayer / 2) + 1);
281 const auto& senseLayer = geo.getSenseLayer(iSLayer);
282 const auto& fieldLayerIn = geo.getFieldLayer(iSLayer - 1);
283 const auto& fieldLayerOut = geo.getFieldLayer(iSLayer);
284
285 rmin_sensitive_left = epLayerBwd.getRmax();
286 rmax_sensitive_left = fieldLayerOut.getR();
287 zback_sensitive_left = senseLayer.getZbwd();
288 zfor_sensitive_left = epLayerBwd.getZfwd();
289
290 rmin_sensitive_middle = fieldLayerIn.getR();
291 rmax_sensitive_middle = fieldLayerOut.getR();
292 zback_sensitive_middle = epLayerBwd.getZfwd();
293 zfor_sensitive_middle = epLayerFwd.getZbwd();
294
295 rmin_sensitive_right = epLayerFwd.getRmax();
296 rmax_sensitive_right = fieldLayerOut.getR();
297 zback_sensitive_right = epLayerFwd.getZbwd();
298 zfor_sensitive_right = senseLayer.getZfwd();
299 } else if (iSLayer >= 7 && iSLayer <= 10) {
300 const auto& epLayerBwd = endplate.getEndPlateLayer(1);
301 const auto& epLayerFwd = endplate.getEndPlateLayer(nEPLayer / 2);
302 const auto& senseLayer = geo.getSenseLayer(iSLayer);
303 const auto& fieldLayerIn = geo.getFieldLayer(iSLayer - 1);
304 const auto& fieldLayerOut = geo.getFieldLayer(iSLayer);
305
306 rmin_sensitive_left = epLayerBwd.getRmax();
307 rmax_sensitive_left = fieldLayerOut.getR();
308 zback_sensitive_left = senseLayer.getZbwd();
309 zfor_sensitive_left = epLayerBwd.getZfwd();
310
311 rmin_sensitive_middle = fieldLayerIn.getR();
312 rmax_sensitive_middle = fieldLayerOut.getR();
313 zback_sensitive_middle = epLayerBwd.getZfwd();
314 zfor_sensitive_middle = epLayerFwd.getZbwd();
315
316 rmin_sensitive_right = epLayerFwd.getRmax();
317 rmax_sensitive_right = fieldLayerOut.getR();
318 zback_sensitive_right = epLayerFwd.getZbwd();
319 zfor_sensitive_right = senseLayer.getZfwd();
320 } else if (iSLayer >= 11 && iSLayer < 47) {
321 const auto& epLayerBwd = endplate.getEndPlateLayer(0);
322 const auto& epLayerFwd = endplate.getEndPlateLayer(nEPLayer / 2);
323 const auto& senseLayer = geo.getSenseLayer(iSLayer);
324 const auto& fieldLayerIn = geo.getFieldLayer(iSLayer - 1);
325 const auto& fieldLayerOut = geo.getFieldLayer(iSLayer);
326
327 rmin_sensitive_left = epLayerBwd.getRmax();
328 rmax_sensitive_left = fieldLayerOut.getR();
329 zback_sensitive_left = senseLayer.getZbwd();
330 zfor_sensitive_left = epLayerBwd.getZfwd();
331
332 rmin_sensitive_middle = fieldLayerIn.getR();
333 rmax_sensitive_middle = fieldLayerOut.getR();
334 zback_sensitive_middle = epLayerBwd.getZfwd();
335 zfor_sensitive_middle = epLayerFwd.getZbwd();
336
337 rmin_sensitive_right = epLayerFwd.getRmax();
338 rmax_sensitive_right = fieldLayerOut.getR();
339 zback_sensitive_right = epLayerFwd.getZbwd();
340 zfor_sensitive_right = senseLayer.getZfwd();
341
342 } else if (iSLayer == 47) {
343
344 const auto& epLayerBwdIn = endplate.getEndPlateLayer(0);
345 const auto& epLayerBwdOut = endplate.getEndPlateLayer((nEPLayer / 2) - 1);
346 const auto& epLayerFwdIn = endplate.getEndPlateLayer(nEPLayer / 2);
347 const auto& epLayerFwdOut = endplate.getEndPlateLayer(nEPLayer - 1);
348 const auto& senseLayer = geo.getSenseLayer(iSLayer);
349
350 int iSLayerMinus1 = iSLayer - 1;
351 const auto& fieldLayerIn = geo.getFieldLayer(iSLayerMinus1);
352 rmin_sensitive_left = epLayerBwdIn.getRmax();
353 rmax_sensitive_left = epLayerBwdOut.getRmax();
354 zback_sensitive_left = senseLayer.getZbwd();
355 zfor_sensitive_left = epLayerBwdIn.getZfwd();
356
357 rmin_sensitive_middle = fieldLayerIn.getR();
358 rmax_sensitive_middle = (geo.getOuterWall(0)).getRmin();
359 zback_sensitive_middle = epLayerBwdIn.getZfwd();
360 zfor_sensitive_middle = epLayerFwdIn.getZbwd();
361
362 rmin_sensitive_right = epLayerFwdIn.getRmax();
363 rmax_sensitive_right = epLayerFwdOut.getRmax();
364 zback_sensitive_right = epLayerFwdIn.getZbwd();
365 zfor_sensitive_right = senseLayer.getZfwd();
366
367 } else {
368 B2ERROR("Undefined sensitive layer : " << iSLayer);
369 continue;
370 }
371
372
373
374 if ((zfor_sensitive_left - zback_sensitive_left) > length_feedthrough) {
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391 G4Tubs* leftTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_leftTube") % iSLayer).str().c_str(),
392 rmin_sensitive_left * CLHEP::cm,
393 rmax_sensitive_left * CLHEP::cm, length_feedthrough * CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
394 G4LogicalVolume* leftTube = new G4LogicalVolume(leftTubeShape, cdcMed,
395 (boost::format("logicalCDCLayer_%1%_leftTube") % iSLayer).str().c_str(), 0, 0, 0);
396 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zback_sensitive_left + length_feedthrough / 2.0)*CLHEP::cm), leftTube,
397 (boost::format(
"physicalCDCLayer_%1%_leftTube") % iSLayer).str().c_str(),
m_logicalCDC,
false, iSLayer);
398
399 G4Tubs* leftSensitiveTubeShape = new G4Tubs((boost::format("solidSD_CDCLayer_%1%_left") % iSLayer).str().c_str(),
400 rmin_sensitive_left * CLHEP::cm, rmax_sensitive_left * CLHEP::cm,
401 (zfor_sensitive_left - zback_sensitive_left - length_feedthrough)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
402 G4LogicalVolume* leftSensitiveTube = new G4LogicalVolume(leftSensitiveTubeShape, cdcMed,
403 (boost::format("logicalSD_CDCLayer_%1%_left") % iSLayer).str().c_str(), 0, 0, 0);
404 leftSensitiveTube->SetSensitiveDetector(
m_sensitive);
405 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_left + zback_sensitive_left + length_feedthrough)*CLHEP::cm / 2.0),
406 leftSensitiveTube, (boost::format(
"physicalSD_CDCLayer_%1%_left") % iSLayer).str().c_str(),
m_logicalCDC,
false, iSLayer);
407 } else {
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424 G4Tubs* leftTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_leftTube") % iSLayer).str().c_str(),
425 rmin_sensitive_left * CLHEP::cm,
426 rmax_sensitive_left * CLHEP::cm, (zfor_sensitive_left - zback_sensitive_left)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
427 G4LogicalVolume* leftTube = new G4LogicalVolume(leftTubeShape, cdcMed,
428 (boost::format("logicalCDCLayer_%1%_leftTube") % iSLayer).str().c_str(), 0, 0, 0);
429 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_left + zback_sensitive_left)*CLHEP::cm / 2.0), leftTube,
430 (boost::format(
"physicalCDCLayer_%1%_leftTube") % iSLayer).str().c_str(),
m_logicalCDC,
false, iSLayer);
431
432
433
434 G4Tubs* leftMidTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_leftMidTube") % iSLayer).str().c_str(),
435 rmin_sensitive_middle * CLHEP::cm, rmax_sensitive_middle * CLHEP::cm,
436 (length_feedthrough - zfor_sensitive_left + zback_sensitive_left)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
437 G4LogicalVolume* leftMidTube = new G4LogicalVolume(leftMidTubeShape, cdcMed,
438 (boost::format("logicalCDCLayer_%1%_leftMidTube") % iSLayer).str().c_str(), 0, 0, 0);
439
440 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (length_feedthrough + zfor_sensitive_left + zback_sensitive_left)*CLHEP::cm / 2.0),
441 leftMidTube, (boost::format(
"physicalCDCLayer_%1%_leftMidTube") % iSLayer).str().c_str(),
m_logicalCDC,
false, iSLayer);
442
443
444 zback_sensitive_middle = length_feedthrough + zback_sensitive_left;
445 }
446
447
448 if ((zfor_sensitive_right - zback_sensitive_right) > length_feedthrough) {
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465 G4Tubs* rightTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_rightTube") % iSLayer).str().c_str(),
466 rmin_sensitive_right * CLHEP::cm, rmax_sensitive_right * CLHEP::cm, length_feedthrough * CLHEP::cm / 2.0, 0 * CLHEP::deg,
467 360.*CLHEP::deg);
468 G4LogicalVolume* rightTube = new G4LogicalVolume(rightTubeShape, cdcMed,
469 (boost::format("logicalCDCLayer_%1%_rightTube") % iSLayer).str().c_str(), 0, 0, 0);
470
471 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_right - length_feedthrough / 2.0)*CLHEP::cm), rightTube,
472 (boost::format(
"physicalCDCLayer_%1%_rightTube") % iSLayer).str().c_str(),
m_logicalCDC,
false, iSLayer);
473
474
475
476 G4Tubs* rightSensitiveTubeShape = new G4Tubs((boost::format("solidSD_CDCLayer_%1%_right") % iSLayer).str().c_str(),
477 rmin_sensitive_right * CLHEP::cm, rmax_sensitive_right * CLHEP::cm,
478 (zfor_sensitive_right - zback_sensitive_right - length_feedthrough)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
479 G4LogicalVolume* rightSensitiveTube = new G4LogicalVolume(rightSensitiveTubeShape, cdcMed,
480 (boost::format("logicalSD_CDCLayer_%1%_right") % iSLayer).str().c_str(), 0, 0, 0);
481 rightSensitiveTube->SetSensitiveDetector(
m_sensitive);
482
483 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_right + zback_sensitive_right - length_feedthrough)*CLHEP::cm / 2.0),
484 rightSensitiveTube, (boost::format(
"physicalSD_CDCLayer_%1%_right") % iSLayer).str().c_str(),
m_logicalCDC,
false, iSLayer);
485
486 } else {
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503 G4Tubs* rightTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_rightTube") % iSLayer).str().c_str(),
504 rmin_sensitive_right * CLHEP::cm, rmax_sensitive_right * CLHEP::cm, (zfor_sensitive_right - zback_sensitive_right)*CLHEP::cm / 2.0,
505 0 * CLHEP::deg, 360.*CLHEP::deg);
506 G4LogicalVolume* rightTube = new G4LogicalVolume(rightTubeShape, cdcMed,
507 (boost::format("logicalCDCLayer_%1%_rightTube") % iSLayer).str().c_str(), 0, 0, 0);
508
509 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_right + zback_sensitive_right)*CLHEP::cm / 2.0), rightTube,
510 (boost::format(
"physicalCDCLayer_%1%_rightTube") % iSLayer).str().c_str(),
m_logicalCDC,
false, iSLayer);
511
512
513
514 G4Tubs* rightMidTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_rightMidTube") % iSLayer).str().c_str(),
515 rmin_sensitive_middle * CLHEP::cm, rmax_sensitive_middle * CLHEP::cm,
516 (length_feedthrough - zfor_sensitive_right + zback_sensitive_right)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
517 G4LogicalVolume* rightMidTube = new G4LogicalVolume(rightMidTubeShape, cdcMed,
518 (boost::format("logicalCDCLayer_%1%_rightMidTube") % iSLayer).str().c_str(), 0, 0, 0);
519 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zback_sensitive_right - length_feedthrough + zfor_sensitive_right)*CLHEP::cm / 2.0),
520 rightMidTube, (boost::format(
"physicalCDCLayer_%1%_rightMidTube") % iSLayer).str().c_str(),
m_logicalCDC,
false, iSLayer);
521
522
523 zfor_sensitive_middle = zfor_sensitive_right - length_feedthrough;
524 }
525
526
527
528 G4Tubs* middleSensitiveTubeShape = new G4Tubs((boost::format("solidSD_CDCLayer_%1%_middle") % iSLayer).str().c_str(),
529 rmin_sensitive_middle * CLHEP::cm, rmax_sensitive_middle * CLHEP::cm,
530 (zfor_sensitive_middle - zback_sensitive_middle)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
531 G4LogicalVolume* middleSensitiveTube = new G4LogicalVolume(middleSensitiveTubeShape, cdcMedGas,
532 (boost::format("logicalSD_CDCLayer_%1%_middle") % iSLayer).str().c_str(), 0, 0, 0);
533
534
535 G4UserLimits* uLimits = new G4UserLimits(8.5 * CLHEP::cm);
537 middleSensitiveTube->SetUserLimits(uLimits);
538 middleSensitiveTube->SetSensitiveDetector(
m_sensitive);
539
540 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_middle + zback_sensitive_middle)*CLHEP::cm / 2.0), middleSensitiveTube,
541 (boost::format(
"physicalSD_CDCLayer_%1%_middle") % iSLayer).str().c_str(),
m_logicalCDC,
false, iSLayer);
542
543
544 if (gcp.getMaterialDefinitionMode() == 2) {
545 G4String sName = "sWire";
546 const G4int jc = 0;
547 B2Vector3D wb0 = cdcgp.wireBackwardPosition(iSLayer, jc);
548
549 B2Vector3D wf0 = cdcgp.wireForwardPosition(iSLayer, jc);
550 G4double tAtZ0 = -wb0.
Z() / (wf0.Z() - wb0.Z());
552
553 const G4double epsl = 25.e-4;
554 G4double reductionBwd = (zback_sensitive_middle + epsl) / wb0.
Z();
555
556 wb0 = reductionBwd * (wb0 - wAtZ0) + wAtZ0;
557
558 G4double reductionFwd = (zfor_sensitive_middle - epsl) / wf0.Z();
559 wf0 = reductionFwd * (wf0 - wAtZ0) + wAtZ0;
560
561 const G4double wireHalfLength = 0.5 * (wf0 - wb0).Mag() * CLHEP::cm;
562 const G4double sWireRadius = 0.5 * cdcgp.senseWireDiameter() * CLHEP::cm;
563
564 G4Tubs* middleSensitiveSwireShape = new G4Tubs(sName, 0., sWireRadius, wireHalfLength, 0., 360. * CLHEP::deg);
565 G4LogicalVolume* middleSensitiveSwire = new G4LogicalVolume(middleSensitiveSwireShape, medTungsten, sName);
566
568
569 G4String fName = "fWire";
570 const G4double fWireRadius = 0.5 * cdcgp.fieldWireDiameter() * CLHEP::cm;
571 G4Tubs* middleSensitiveFwireShape = new G4Tubs(fName, 0., fWireRadius, wireHalfLength, 0., 360. * CLHEP::deg);
572 G4LogicalVolume* middleSensitiveFwire = new G4LogicalVolume(middleSensitiveFwireShape, medAluminum, fName);
573
575
576 const G4double diameter = cdcgp.fieldWireDiameter();
577
578 const G4int nCells = cdcgp.nWiresInLayer(iSLayer);
579 const G4double dphi = M_PI / nCells;
581
582 for (int ic = 0; ic < nCells; ++ic) {
583
584 B2Vector3D wb = cdcgp.wireBackwardPosition(iSLayer, ic);
585 B2Vector3D wf = cdcgp.wireForwardPosition(iSLayer, ic);
586 G4double tAtZ02 = -wb.Z() / (wf.Z() - wb.Z());
588 G4double reductionBwd2 = (zback_sensitive_middle + epsl) / wb.Z();
589 wb = reductionBwd2 * (wb - wAtZ02) + wAtZ02;
590 G4double reductionFwd2 = (zfor_sensitive_middle - epsl) / wf.Z();
591 wf = reductionFwd2 * (wf - wAtZ02) + wAtZ02;
592
593 G4double thetaYZ = -asin((wf - wb).Y() / (wf - wb).Mag());
594
595 B2Vector3D fMinusBInZX((wf - wb).X(), 0., (wf - wb).Z());
596 G4double thetaZX = asin((unitZ.Cross(fMinusBInZX)).Y() / fMinusBInZX.Mag());
597 G4RotationMatrix rotM;
598
599 rotM.rotateX(thetaYZ * CLHEP::rad);
600 rotM.rotateY(thetaZX * CLHEP::rad);
601
602 G4ThreeVector xyz(0.5 * (wb.X() + wf.X()) * CLHEP::cm,
603 0.5 * (wb.Y() + wf.Y()) * CLHEP::cm, 0.);
604
605
606
607 new G4PVPlacement(G4Transform3D(rotM, xyz), middleSensitiveSwire, sName, middleSensitiveTube, false, ic);
608
609
611 G4double rF = rmax_sensitive_middle - 0.5 * diameter;
612
613 G4double phi = atan2(wbF.Y(), wbF.X());
614 wbF.SetX(rF * cos(phi));
615 wbF.SetY(rF * sin(phi));
616
618 rF = rmax_sensitive_middle - 0.5 * diameter;
619 phi = atan2(wfF.Y(), wfF.X());
620 wfF.SetX(rF * cos(phi));
621 wfF.SetY(rF * sin(phi));
622
623 thetaYZ = -asin((wfF - wbF).Y() / (wfF - wbF).Mag());
624
625 fMinusBInZX = wfF - wbF;
626 fMinusBInZX.SetY(0.);
627 thetaZX = asin((unitZ.Cross(fMinusBInZX)).Y() / fMinusBInZX.Mag());
628
629 G4RotationMatrix rotM1;
630 rotM1.rotateX(thetaYZ * CLHEP::rad);
631 rotM1.rotateY(thetaZX * CLHEP::rad);
632
633 xyz.setX(0.5 * (wbF.X() + wfF.X()) * CLHEP::cm);
634 xyz.setY(0.5 * (wbF.Y() + wfF.Y()) * CLHEP::cm);
635
636 if (iSLayer != nSLayer - 1) {
637
638 new G4PVPlacement(G4Transform3D(rotM1, xyz), middleSensitiveFwire, fName, middleSensitiveTube, false, ic);
639 }
640
641
642 wbF = wb;
643 rF = wbF.Perp();
644 phi = atan2(wbF.Y(), wbF.X());
645 wbF.SetX(rF * cos(phi + dphi));
646 wbF.SetY(rF * sin(phi + dphi));
647
648 wfF = wf;
649 rF = wfF.Perp();
650 phi = atan2(wfF.Y(), wfF.X());
651 wfF.SetX(rF * cos(phi + dphi));
652 wfF.SetY(rF * sin(phi + dphi));
653
654 thetaYZ = -asin((wfF - wbF).Y() / (wfF - wbF).Mag());
655
656 fMinusBInZX = wfF - wbF;
657 fMinusBInZX.SetY(0.);
658 thetaZX = asin((unitZ.Cross(fMinusBInZX)).Y() / fMinusBInZX.Mag());
659
660 G4RotationMatrix rotM2;
661 rotM2.rotateX(thetaYZ * CLHEP::rad);
662 rotM2.rotateY(thetaZX * CLHEP::rad);
663
664 xyz.setX(0.5 * (wbF.X() + wfF.X()) * CLHEP::cm);
665 xyz.setY(0.5 * (wbF.Y() + wfF.Y()) * CLHEP::cm);
666
667
668 new G4PVPlacement(G4Transform3D(rotM2, xyz), middleSensitiveFwire, fName, middleSensitiveTube, false, ic + nCells);
669
670
671 wbF = wb;
672 rF = rmax_sensitive_middle - 0.5 * diameter;
673 phi = atan2(wbF.Y(), wbF.X());
674 wbF.SetX(rF * cos(phi + dphi));
675 wbF.SetY(rF * sin(phi + dphi));
676
677 wfF = wf;
678 rF = rmax_sensitive_middle - 0.5 * diameter;
679 phi = atan2(wfF.Y(), wfF.X());
680 wfF.SetX(rF * cos(phi + dphi));
681 wfF.SetY(rF * sin(phi + dphi));
682
683 thetaYZ = -asin((wfF - wbF).Y() / (wfF - wbF).Mag());
684
685 fMinusBInZX = wfF - wbF;
686 fMinusBInZX.SetY(0.);
687 thetaZX = asin((unitZ.Cross(fMinusBInZX)).Y() / fMinusBInZX.Mag());
688
689 G4RotationMatrix rotM3;
690 rotM3.rotateX(thetaYZ * CLHEP::rad);
691 rotM3.rotateY(thetaZX * CLHEP::rad);
692
693 xyz.setX(0.5 * (wbF.X() + wfF.X()) * CLHEP::cm);
694 xyz.setY(0.5 * (wbF.Y() + wfF.Y()) * CLHEP::cm);
695
696 if (iSLayer != nSLayer - 1) {
697 new G4PVPlacement(G4Transform3D(rotM3, xyz), middleSensitiveFwire, fName, middleSensitiveTube, false, ic + 2 * nCells);
698 }
699 }
700 }
701
702 }
703
704
705
706
707 m_VisAttributes.push_back(
new G4VisAttributes(
true, G4Colour(1., 1., 0.)));
708 for (const auto& endplate : geo.getEndPlates()) {
709 for (const auto& epLayer : endplate.getEndPlateLayers()) {
710 const int iEPLayer = epLayer.getILayer();
711 const string name = epLayer.getName();
712 const double rmin = epLayer.getRmin();
713 const double rmax = epLayer.getRmax();
714 const double zbwd = epLayer.getZbwd();
715 const double zfwd = epLayer.getZfwd();
716 const double length = (zfwd - zbwd) / 2.0;
717
718 G4Tubs* tube = new G4Tubs("solidCDCEndplate" + name, rmin * CLHEP::cm,
719 rmax * CLHEP::cm, length * CLHEP::cm, 0 * CLHEP::deg, 360.*CLHEP::deg);
720 G4LogicalVolume* logical =
new G4LogicalVolume(tube,
Materials::get(
"G4_Al"),
721 "logicalCDCEndplate" + name, 0, 0);
723 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfwd + zbwd)*CLHEP::cm / 2.0), logical,
724 "physicalCDCEndplate" + name,
m_logicalCDC,
false, iEPLayer);
725
726 }
727 }
728
729
730
731 for (const auto& frontend : geo.getFrontends()) {
732
733 const int iEB = frontend.getId();
734 const double ebInnerR = frontend.getRmin();
735 const double ebOuterR = frontend.getRmax();
736 const double ebBZ = frontend.getZbwd();
737 const double ebFZ = frontend.getZfwd();
738
739 G4Tubs* ebTubeShape = new G4Tubs((boost::format("solidSD_ElectronicsBoard_Layer%1%") % iEB).str().c_str(), ebInnerR * CLHEP::cm,
740 ebOuterR * CLHEP::cm, (ebFZ - ebBZ)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
741
742 G4LogicalVolume* ebTube = new G4LogicalVolume(ebTubeShape, medNEMA_G10_Plate,
743 (boost::format("logicalSD_ElectronicsBoard_Layer%1%") % iEB).str().c_str(), 0, 0, 0);
747 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (ebFZ + ebBZ)*CLHEP::cm / 2.0), ebTube,
748 (boost::format(
"physicalSD_ElectronicsBoard_Layer%1%") % iEB).str().c_str(),
m_logicalCDC,
false, iEB);
749 }
750
751
752
753
755
756
757
758
760
761
762
763
765
766
767
768
769 for (const auto& rib : geo.getRibs()) {
770
771 const int id = rib.getId();
772 const double length = rib.getLength();
773 const double width = rib.getWidth();
774 const double thick = rib.getThick();
775 const double rotx = rib.getRotX();
776 const double roty = rib.getRotY();
777 const double rotz = rib.getRotZ();
778 const double x = rib.getX();
779 const double y = rib.getY();
780 const double z = rib.getZ();
781 const int offset = rib.getOffset();
782 const int ndiv = rib.getNDiv();
783
784 const string solidName = "solidRib" + to_string(id);
785 const string logicalName = "logicalRib" + to_string(id);
786 G4Box* boxShape = new G4Box(solidName, 0.5 * length * CLHEP::cm,
787 0.5 * width * CLHEP::cm,
788 0.5 * thick * CLHEP::cm);
789
790 const double rmax = 0.5 * length;
791 const double rmin = max((rmax - thick), 0.);
792 G4Tubs* tubeShape = new G4Tubs(solidName,
793 rmin * CLHEP::cm,
794 rmax * CLHEP::cm,
795 0.5 * width * CLHEP::cm,
796 0.,
797 360. * CLHEP::deg);
798
799
800
801
802 G4LogicalVolume* logicalV = new G4LogicalVolume(boxShape, medAluminum, logicalName, 0, 0, 0);
803 if (id > 39 && id < 78)
804 logicalV = new G4LogicalVolume(boxShape, medCopper, logicalName, 0, 0, 0);
805 if ((id > 77 && id < 94) || (id > 131 && id < 146))
806 logicalV = new G4LogicalVolume(boxShape, medNEMA_G10_Plate, logicalName, 0, 0, 0);
807 if (id > 93 && id < 110)
808 logicalV = new G4LogicalVolume(tubeShape, medCopper, logicalName, 0, 0, 0);
809 if (id > 109 && id < 126)
810 logicalV = new G4LogicalVolume(tubeShape, medH2O, logicalName, 0, 0, 0);
811 if (id > 127 && id < 132)
812 logicalV = new G4LogicalVolume(boxShape, medHV, logicalName, 0, 0, 0);
813
814
815
816
817
818
819
821
822 const double phi = 360.0 / ndiv;
823
824 G4RotationMatrix rot = G4RotationMatrix();
825 double dz = thick;
826 if (id > 93 && id < 126) dz = 0;
827
828 G4ThreeVector arm(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - dz * CLHEP::cm / 2.0);
829 rot.rotateX(rotx);
830 rot.rotateY(roty);
831 rot.rotateZ(rotz);
832 if (offset) {
833 rot.rotateZ(0.5 * phi * CLHEP::deg);
834 arm.rotateZ(0.5 * phi * CLHEP::deg);
835 }
836 for (int i = 0; i < ndiv; ++i) {
837 const string physicalName = "physicalRib_" + to_string(id) + " " + to_string(i);
838 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
840 rot.rotateZ(phi * CLHEP::deg);
841 arm.rotateZ(phi * CLHEP::deg);
842 }
843
844 }
845
846
847
848
849 for (const auto& rib2 : geo.getRib2s()) {
850
851 const int id = rib2.getId();
852 const double length = rib2.getLength();
853 const double width = rib2.getWidth();
854 const double thick = rib2.getThick();
855 const double width2 = rib2.getWidth2();
856 const double thick2 = rib2.getThick2();
857 const double rotx = rib2.getRotX();
858 const double roty = rib2.getRotY();
859 const double rotz = rib2.getRotZ();
860 const double x = rib2.getX();
861 const double y = rib2.getY();
862 const double z = rib2.getZ();
863 const int ndiv = rib2.getNDiv();
864
865 const string solidName = "solidRib2" + to_string(id);
866 const string logicalName = "logicalRib2" + to_string(id);
867 G4Trd* trdShape = new G4Trd(solidName,
868 0.5 * thick * CLHEP::cm,
869 0.5 * thick2 * CLHEP::cm,
870 0.5 * width * CLHEP::cm,
871 0.5 * width2 * CLHEP::cm,
872 0.5 * length * CLHEP::cm);
873
874 G4LogicalVolume* logicalV = new G4LogicalVolume(trdShape, medAluminum, logicalName, 0, 0, 0);
875
876 if (id > 0)
877 logicalV = new G4LogicalVolume(trdShape, medCopper, logicalName, 0, 0, 0);
878
880
881 const double phi = 360.0 / ndiv;
882
883 G4RotationMatrix rot = G4RotationMatrix();
884 G4ThreeVector arm(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - thick * CLHEP::cm / 2.0);
885
886 rot.rotateX(rotx);
887 rot.rotateY(roty);
888 rot.rotateZ(rotz);
889 for (int i = 0; i < ndiv; ++i) {
890 const string physicalName = "physicalRib2_" + to_string(id) + " " + to_string(i);
891 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
893 rot.rotateZ(phi * CLHEP::deg);
894 arm.rotateZ(phi * CLHEP::deg);
895 }
896
897 }
898
899
900
901
902 for (const auto& rib3 : geo.getRib3s()) {
903
904 const int id = rib3.getId();
905 const double length = rib3.getLength();
906 const double width = rib3.getWidth();
907 const double thick = rib3.getThick();
908 const double r = rib3.getR();
909 const double x = rib3.getX();
910 const double y = rib3.getY();
911 const double z = rib3.getZ();
912 const double rx = rib3.getRx();
913 const double ry = rib3.getRy();
914 const double rz = rib3.getRz();
915 const int offset = rib3.getOffset();
916 const int ndiv = rib3.getNDiv();
917
918 const string logicalName = "logicalRib3" + to_string(id);
919 G4VSolid* boxShape = new G4Box("Block",
920 0.5 * length * CLHEP::cm,
921 0.5 * width * CLHEP::cm,
922 0.5 * thick * CLHEP::cm);
923 G4VSolid* tubeShape = new G4Tubs("Hole",
924 0.,
925 r * CLHEP::cm,
926 width * CLHEP::cm,
927 0. * CLHEP::deg,
928 360. * CLHEP::deg);
929
930 G4RotationMatrix rotsub = G4RotationMatrix();
931 rotsub.rotateX(90. * CLHEP::deg);
932 G4ThreeVector trnsub(rx * CLHEP::cm - x * CLHEP::cm, ry * CLHEP::cm - y * CLHEP::cm,
933 rz * CLHEP::cm - z * CLHEP::cm + 0.5 * thick * CLHEP::cm);
934 G4VSolid* coolingBlock = new G4SubtractionSolid("Block-Hole",
935 boxShape,
936 tubeShape,
937 G4Transform3D(rotsub,
938 trnsub));
939
940 G4LogicalVolume* logicalV = new G4LogicalVolume(coolingBlock, medCopper, logicalName, 0, 0, 0);
941
943
944 const double phi = 360.0 / ndiv;
945
946 G4RotationMatrix rot = G4RotationMatrix();
947 G4ThreeVector arm(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - thick * CLHEP::cm / 2.0);
948
949 if (offset) {
950 rot.rotateZ(0.5 * phi * CLHEP::deg);
951 arm.rotateZ(0.5 * phi * CLHEP::deg);
952 }
953 for (int i = 0; i < ndiv; ++i) {
954 const string physicalName = "physicalRib3_" + to_string(id) + " " + to_string(i);
955 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
957 rot.rotateZ(phi * CLHEP::deg);
958 arm.rotateZ(phi * CLHEP::deg);
959 }
960
961 }
962
963
964
965
966 for (const auto& rib4 : geo.getRib4s()) {
967
968 const int id = rib4.getId();
969 const double length = rib4.getLength();
970 const double width = rib4.getWidth();
971 const double thick = rib4.getThick();
972 const double length2 = rib4.getLength2();
973 const double width2 = rib4.getWidth2();
974 const double thick2 = rib4.getThick2();
975 const double x = rib4.getX();
976 const double y = rib4.getY();
977 const double z = rib4.getZ();
978 const double x2 = rib4.getX2();
979 const double y2 = rib4.getY2();
980 const double z2 = rib4.getZ2();
981 const int offset = rib4.getOffset();
982 const int ndiv = rib4.getNDiv();
983
984 const string logicalName = "logicalRib4" + to_string(id);
985 G4VSolid* baseShape = new G4Box("Base",
986 0.5 * length * CLHEP::cm,
987 0.5 * width * CLHEP::cm,
988 0.5 * thick * CLHEP::cm);
989 G4VSolid* sqShape = new G4Box("Sq",
990 0.5 * length2 * CLHEP::cm,
991 0.5 * width2 * CLHEP::cm,
992 0.5 * thick2 * CLHEP::cm);
993
994 G4RotationMatrix rotsub = G4RotationMatrix();
995 double dzc = (z2 - thick2 / 2.) - (z - thick / 2.);
996 G4ThreeVector trnsub(x2 * CLHEP::cm - x * CLHEP::cm,
997 y2 * CLHEP::cm - y * CLHEP::cm,
998 dzc * CLHEP::cm);
999 G4VSolid* sqHoleBase = new G4SubtractionSolid("Box-Sq",
1000 baseShape,
1001 sqShape,
1002 G4Transform3D(rotsub,
1003 trnsub)
1004 );
1005
1006 G4LogicalVolume* logicalV = new G4LogicalVolume(sqHoleBase, medCopper, logicalName, 0, 0, 0);
1007 if (id < 19) {
1008 logicalV = new G4LogicalVolume(sqHoleBase, medNEMA_G10_Plate, logicalName, 0, 0, 0);
1009 logicalV->SetSensitiveDetector(new BkgSensitiveDetector("CDC", 2000 + id));
1010 }
1011
1013
1014 const double phi = 360.0 / ndiv;
1015
1016 G4RotationMatrix rot = G4RotationMatrix();
1017 G4ThreeVector arm(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - thick * CLHEP::cm / 2.0);
1018
1019 if (offset) {
1020 rot.rotateZ(0.5 * phi * CLHEP::deg);
1021 arm.rotateZ(0.5 * phi * CLHEP::deg);
1022 }
1023 for (int i = 0; i < ndiv; ++i) {
1024 const string physicalName = "physicalRib4_" + to_string(id) + " " + to_string(i);
1025 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
1027 rot.rotateZ(phi * CLHEP::deg);
1028 arm.rotateZ(phi * CLHEP::deg);
1029 }
1030
1031 }
1032
1033
1034
1035 for (const auto& rib5 : geo.getRib5s()) {
1036
1037 const int id = rib5.getId();
1038 const double dr = rib5.getDr();
1039 const double dz = rib5.getDz();
1040 const double width = rib5.getWidth();
1041 const double thick = rib5.getThick();
1042 const double rin = rib5.getRin();
1043 const double x = rib5.getX();
1044 const double y = rib5.getY();
1045 const double z = rib5.getZ();
1046 const double rotx = rib5.getRotx();
1047 const double roty = rib5.getRoty();
1048 const double rotz = rib5.getRotz();
1049 const int offset = rib5.getOffset();
1050 const int ndiv = rib5.getNDiv();
1051
1052 const string solidName = "solidRib5" + to_string(id);
1053 const string logicalName = "logicalRib5" + to_string(id);
1054
1055 const double rmax = rin + thick;
1056 const double rmin = rin;
1057 const double dphi = 2. * atan2(dz, dr);
1058 const double ddphi = thick * tan(dphi) / rin;
1059 const double ddphi2 = width / 2. * width / 2. / (x + dr) / rin;
1060 const double cphi = dphi - ddphi - ddphi2;
1061 G4Tubs* tubeShape = new G4Tubs(solidName,
1062 rmin * CLHEP::cm,
1063 rmax * CLHEP::cm,
1064 0.5 * width * CLHEP::cm,
1065 0.,
1066 cphi);
1067
1068 G4LogicalVolume* logicalV = new G4LogicalVolume(tubeShape, medAluminum, logicalName, 0, 0, 0);
1069
1071
1072 const double phi = 360.0 / ndiv;
1073
1074 G4RotationMatrix rot = G4RotationMatrix();
1075
1076
1077 G4ThreeVector arm(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - rin * CLHEP::cm - thick * CLHEP::cm);
1078 rot.rotateX(rotx);
1079 rot.rotateY(roty);
1080 rot.rotateZ(rotz);
1081 if (offset) {
1082 rot.rotateZ(0.5 * phi * CLHEP::deg);
1083 arm.rotateZ(0.5 * phi * CLHEP::deg);
1084 }
1085 for (int i = 0; i < ndiv; ++i) {
1086 const string physicalName = "physicalRib5_" + to_string(id) + " " + to_string(i);
1087 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
1089 rot.rotateZ(phi * CLHEP::deg);
1090 arm.rotateZ(phi * CLHEP::deg);
1091 }
1092
1093 }
1094
1095
1097 }
DataType Z() const
access variable Z (= .at(2) without boundary check)
static CDCGeometryPar & Instance(const CDCGeometry *=nullptr)
Static method to get a reference to the CDCGeometryPar instance.
void createNeutronShields(const GearDir &content)
Create neutron shield from gearbox.
void createCover2s(const GearDir &content)
Create CDC cover2s from gear box.
void createMapper(G4LogicalVolume &topVolume)
Create the B-field mapper geometry (tentative function)
void createCovers(const GearDir &content)
Create CDC covers from gear box.
B2Vector3< double > B2Vector3D
typedef for common usage with double