Belle II Software development
GeoCDCCreator.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#include <cdc/geometry/GeoCDCCreator.h>
10#include <cdc/geometry/CDCGeometryPar.h>
11#include <cdc/geometry/CDCGeoControlPar.h>
12#include <cdc/simulation/CDCSimControlPar.h>
13#include <cdc/simulation/CDCSensitiveDetector.h>
14#include <simulation/background/BkgSensitiveDetector.h>
15
16#include <geometry/CreatorFactory.h>
17#include <geometry/Materials.h>
18
19#include <cmath>
20#include <boost/format.hpp>
21
22#include <G4Material.hh>
23#include <G4Box.hh>
24#include <G4Tubs.hh>
25#include <G4Torus.hh>
26#include <G4Trd.hh>
27#include <G4SubtractionSolid.hh>
28#include <G4Region.hh>
29#include <G4VSolid.hh>
30
31#include <G4Polycone.hh>
32#include <G4Cons.hh>
33#include <G4Colour.hh>
34#include <G4LogicalVolume.hh>
35#include <G4PVPlacement.hh>
36#include <G4Transform3D.hh>
37#include <G4VisAttributes.hh>
38#include <G4RotationMatrix.hh>
39#include <G4UserLimits.hh>
40#include <iostream>
41
42using namespace std;
43
44namespace Belle2 {
49
50 using namespace geometry;
51
52 namespace CDC {
53
57
59
60 //-----------------------------------------------------------------
61 // Implementation
62 //-----------------------------------------------------------------
63
65 {
66 // Set job control params. before sensitivedetector and geometry construction
69
70 m_logicalCDC = 0;
71 m_physicalCDC = 0;
72 m_VisAttributes.clear();
73 m_VisAttributes.push_back(new G4VisAttributes(false)); // for "invisible"
74 m_userLimits.clear();
75 }
76
77
79 {
80 delete m_sensitive;
82 for (BkgSensitiveDetector* sensitiveDetector : m_BkgSensitiveRib4)
83 delete sensitiveDetector;
84 for (G4VisAttributes* visAttr : m_VisAttributes) delete visAttr;
85 m_VisAttributes.clear();
86 for (G4UserLimits* userLimits : m_userLimits) delete userLimits;
87 m_userLimits.clear();
88 }
89
90 void GeoCDCCreator::createGeometry(const CDCGeometry& geo, G4LogicalVolume& topVolume, geometry::GeometryTypes)
91 {
92
93 m_sensitive = new CDCSensitiveDetector("CDCSensitiveDetector", (2 * 24)* CLHEP::eV, 10 * CLHEP::MeV);
94
95 const G4double realTemperture = (273.15 + 23.) * CLHEP::kelvin;
96 G4Material* medHelium = geometry::Materials::get("CDCHeGas");
97 G4Material* medEthane = geometry::Materials::get("CDCEthaneGas");
98 G4Material* medAluminum = geometry::Materials::get("Al");
99 G4Material* medTungsten = geometry::Materials::get("W");
100 G4Material* medCFRP = geometry::Materials::get("CFRP");
101 G4Material* medNEMA_G10_Plate = geometry::Materials::get("NEMA_G10_Plate");
102 G4Material* medGlue = geometry::Materials::get("CDCGlue");
103 G4Material* medAir = geometry::Materials::get("Air");
104
105 G4double h2odensity = 1.000 * CLHEP::g / CLHEP::cm3;
106 G4double a = 1.01 * CLHEP::g / CLHEP::mole;
107 G4Element* elH = new G4Element("Hydrogen", "H", 1., a);
108 a = 16.00 * CLHEP::g / CLHEP::mole;
109 G4Element* elO = new G4Element("Oxygen", "O", 8., a);
110 G4Material* medH2O = new G4Material("Water", h2odensity, 2);
111 medH2O->AddElement(elH, 2);
112 medH2O->AddElement(elO, 1);
113 G4Material* medCopper = geometry::Materials::get("Cu");
114 G4Material* medHV = geometry::Materials::get("CDCHVCable");
115 //G4Material* medFiber = geometry::Materials::get("CDCOpticalFiber");
116 //G4Material* medCAT7 = geometry::Materials::get("CDCCAT7");
117 //G4Material* medTRG = geometry::Materials::get("CDCOpticalFiberTRG");
118
119 // Total cross section
120 const double rmax_innerWall = geo.getFiducialRmin();
121 const double rmin_outerWall = geo.getFiducialRmax();
122 const double diameter_senseWire = geo.getSenseDiameter();
123 const double diameter_fieldWire = geo.getFieldDiameter();
124 const double num_senseWire = static_cast<double>(geo.getNSenseWires());
125 const double num_fieldWire = static_cast<double>(geo.getNFieldWires());
126 double totalCS = M_PI * (rmin_outerWall * rmin_outerWall - rmax_innerWall * rmax_innerWall);
127
128 // Sense wire cross section
129 double senseCS = M_PI * (diameter_senseWire / 2) * (diameter_senseWire / 2) * num_senseWire;
130
131 // Field wire cross section
132 double fieldCS = M_PI * (diameter_fieldWire / 2) * (diameter_fieldWire / 2) * num_fieldWire;
133
134 // Density
135 const double denHelium = medHelium->GetDensity() / 2.0;
136 const double denEthane = medEthane->GetDensity() / 2.0;
137 const double denAluminum = medAluminum->GetDensity() * (fieldCS / totalCS);
138 const double denTungsten = medTungsten->GetDensity() * (senseCS / totalCS);
139 const double density = denHelium + denEthane + denAluminum + denTungsten;
140 G4Material* cdcMed = new G4Material("CDCGasWire", density, 4, kStateGas, realTemperture);
141 cdcMed->AddMaterial(medHelium, denHelium / density);
142 cdcMed->AddMaterial(medEthane, denEthane / density);
143 cdcMed->AddMaterial(medTungsten, denTungsten / density);
144 cdcMed->AddMaterial(medAluminum, denAluminum / density);
145
146 G4Material* cdcMedGas = cdcMed;
147
150 // std::cout << gcp.getMaterialDefinitionMode() << std::endl;
151
152 // if (cdcgp.getMaterialDefinitionMode() == 2) {
153 if (gcp.getMaterialDefinitionMode() == 2) {
154 const double density2 = denHelium + denEthane;
155 cdcMedGas = new G4Material("CDCRealGas", density2, 2, kStateGas, realTemperture);
156 cdcMedGas->AddMaterial(medHelium, denHelium / density2);
157 cdcMedGas->AddMaterial(medEthane, denEthane / density2);
158 }
159
160 if (gcp.getPrintMaterialTable()) {
161 G4cout << *(G4Material::GetMaterialTable());
162 }
163
164 const auto& mother = geo.getMotherVolume();
165 const auto& motherRmin = mother.getRmin();
166 const auto& motherRmax = mother.getRmax();
167 const auto& motherZ = mother.getZ();
168 G4Polycone* solid_cdc =
169 new G4Polycone("solidCDC", 0 * CLHEP::deg, 360.* CLHEP::deg,
170 mother.getNNodes(), motherZ.data(),
171 motherRmin.data(), motherRmax.data());
172 m_logicalCDC = new G4LogicalVolume(solid_cdc, medAir, "logicalCDC", 0, 0, 0);
173 m_physicalCDC = new G4PVPlacement(0, G4ThreeVector(geo.getGlobalOffsetX() * CLHEP::cm,
174 geo.getGlobalOffsetY() * CLHEP::cm,
175 geo.getGlobalOffsetZ() * CLHEP::cm), m_logicalCDC,
176 "physicalCDC", &topVolume, false, 0);
177
178 // Set up region for production cuts
179 G4Region* aRegion = new G4Region("CDCEnvelope");
180 m_logicalCDC->SetRegion(aRegion);
181 aRegion->AddRootLogicalVolume(m_logicalCDC);
182
183 m_VisAttributes.push_back(new G4VisAttributes(true, G4Colour(0., 1., 0.)));
184 for (const auto& wall : geo.getOuterWalls()) {
185 const int iOuterWall = wall.getId();
186 const string wallName = wall.getName();
187 const double wallRmin = wall.getRmin();
188 const double wallRmax = wall.getRmax();
189 const double wallZfwd = wall.getZfwd();
190 const double wallZbwd = wall.getZbwd();
191 const double length = (wallZfwd - wallZbwd) / 2.0;
192
193
194 G4Material* medWall;
195 if (strstr((wallName).c_str(), "MiddleWall") != nullptr) {
196 medWall = medCFRP;
197 } else {
198 medWall = medAluminum;
199 }
200 G4Tubs* outerWallTubeShape = new G4Tubs("solid" + wallName, wallRmin * CLHEP::cm,
201 wallRmax * CLHEP::cm, length * CLHEP::cm, 0 * CLHEP::deg, 360.*CLHEP::deg);
202
203 G4LogicalVolume* outerWallTube = new G4LogicalVolume(outerWallTubeShape, medWall, "solid" + wallName, 0, 0, 0);
204 outerWallTube->SetVisAttributes(m_VisAttributes.back());
205 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (length + wallZbwd)*CLHEP::cm), outerWallTube, "logical" + wallName,
206 m_logicalCDC, false, iOuterWall);
207 }
208
209
210 m_VisAttributes.push_back(new G4VisAttributes(true, G4Colour(0., 1., 0.)));
211 for (const auto& wall : geo.getInnerWalls()) {
212 const string wallName = wall.getName();
213 const double wallRmin = wall.getRmin();
214 const double wallRmax = wall.getRmax();
215 const double wallZfwd = wall.getZfwd();
216 const double wallZbwd = wall.getZbwd();
217 const double length = (wallZfwd - wallZbwd) / 2.0;
218 const int iInnerWall = wall.getId();
219
220 G4Material* medWall;
221 if (strstr(wallName.c_str(), "MiddleWall") != nullptr) {
222 medWall = medCFRP;
223 } else if (strstr(wallName.c_str(), "MiddleGlue") != nullptr) { // Glue layer 0.005 mmt
224 medWall = medGlue;
225 } else { // Al layer 0.1 mmt
226 medWall = medAluminum;
227 }
228
229 G4Tubs* innerWallTubeShape = new G4Tubs("solid" + wallName, wallRmin * CLHEP::cm,
230 wallRmax * CLHEP::cm, length * CLHEP::cm, 0 * CLHEP::deg, 360.*CLHEP::deg);
231 G4LogicalVolume* innerWallTube = new G4LogicalVolume(innerWallTubeShape, medWall, "logical" + wallName, 0, 0, 0);
232 innerWallTube->SetVisAttributes(m_VisAttributes.back());
233 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (length + wallZbwd)*CLHEP::cm), innerWallTube, "physical" + wallName,
234 m_logicalCDC, false, iInnerWall);
235
236
237 }
238
239
240
241 //
242 // Construct sensitive layers.
243 //
244 const int nSLayer = geo.getNSenseLayers();
245 const double length_feedthrough = geo.getFeedthroughLength();
246 for (int iSLayer = 0; iSLayer < nSLayer; ++iSLayer) {
247 const auto& endplate = geo.getEndPlate(iSLayer);
248 const int nEPLayer = endplate.getNEndPlateLayers();
249 // Get parameters for sensitive layer: left, middle and right.
250 double rmin_sensitive_left, rmax_sensitive_left;
251 double rmin_sensitive_middle, rmax_sensitive_middle;
252 double rmin_sensitive_right, rmax_sensitive_right;
253 double zback_sensitive_left, zfor_sensitive_left;
254 double zback_sensitive_middle, zfor_sensitive_middle;
255 double zback_sensitive_right, zfor_sensitive_right;
256
257 if (iSLayer == 0) {
258 const auto& epLayerBwd = endplate.getEndPlateLayer(0);
259 const auto& epLayerFwd = endplate.getEndPlateLayer(nEPLayer / 2);
260 const auto& senseLayer = geo.getSenseLayer(iSLayer);
261 const auto& fieldLayer = geo.getFieldLayer(iSLayer);
262
263 rmin_sensitive_left = epLayerBwd.getRmax();
264 rmax_sensitive_left = fieldLayer.getR();
265 zback_sensitive_left = senseLayer.getZbwd();
266 zfor_sensitive_left = epLayerBwd.getZfwd();
267
268 rmin_sensitive_middle = (geo.getInnerWall(0)).getRmax();
269 rmax_sensitive_middle = fieldLayer.getR();
270 zback_sensitive_middle = epLayerBwd.getZfwd();
271 zfor_sensitive_middle = epLayerFwd.getZbwd();
272
273 rmin_sensitive_right = epLayerFwd.getRmax();
274 rmax_sensitive_right = fieldLayer.getR();
275 zback_sensitive_right = epLayerFwd.getZbwd();
276 zfor_sensitive_right = senseLayer.getZfwd();
277 } else if (iSLayer >= 1 && iSLayer <= 14) {
278 const auto& epLayerBwd = endplate.getEndPlateLayer(1);
279 const auto& epLayerFwd = endplate.getEndPlateLayer((nEPLayer / 2) + 1);
280 const auto& senseLayer = geo.getSenseLayer(iSLayer);
281 const auto& fieldLayerIn = geo.getFieldLayer(iSLayer - 1);
282 const auto& fieldLayerOut = geo.getFieldLayer(iSLayer);
283
284 rmin_sensitive_left = epLayerBwd.getRmax();
285 rmax_sensitive_left = fieldLayerOut.getR();
286 zback_sensitive_left = senseLayer.getZbwd();
287 zfor_sensitive_left = epLayerBwd.getZfwd();
288
289 rmin_sensitive_middle = fieldLayerIn.getR();
290 rmax_sensitive_middle = fieldLayerOut.getR();
291 zback_sensitive_middle = epLayerBwd.getZfwd();
292 zfor_sensitive_middle = epLayerFwd.getZbwd();
293
294 rmin_sensitive_right = epLayerFwd.getRmax();
295 rmax_sensitive_right = fieldLayerOut.getR();
296 zback_sensitive_right = epLayerFwd.getZbwd();
297 zfor_sensitive_right = senseLayer.getZfwd();
298 } else if (iSLayer >= 15 && iSLayer <= 18) {
299 const auto& epLayerBwd = endplate.getEndPlateLayer(1);
300 const auto& epLayerFwd = endplate.getEndPlateLayer(nEPLayer / 2);
301 const auto& senseLayer = geo.getSenseLayer(iSLayer);
302 const auto& fieldLayerIn = geo.getFieldLayer(iSLayer - 1);
303 const auto& fieldLayerOut = geo.getFieldLayer(iSLayer);
304
305 rmin_sensitive_left = epLayerBwd.getRmax();
306 rmax_sensitive_left = fieldLayerOut.getR();
307 zback_sensitive_left = senseLayer.getZbwd();
308 zfor_sensitive_left = epLayerBwd.getZfwd();
309
310 rmin_sensitive_middle = fieldLayerIn.getR();
311 rmax_sensitive_middle = fieldLayerOut.getR();
312 zback_sensitive_middle = epLayerBwd.getZfwd();
313 zfor_sensitive_middle = epLayerFwd.getZbwd();
314
315 rmin_sensitive_right = epLayerFwd.getRmax();
316 rmax_sensitive_right = fieldLayerOut.getR();
317 zback_sensitive_right = epLayerFwd.getZbwd();
318 zfor_sensitive_right = senseLayer.getZfwd();
319 } else if (iSLayer >= 19 && iSLayer < 55) {
320 const auto& epLayerBwd = endplate.getEndPlateLayer(0);
321 const auto& epLayerFwd = endplate.getEndPlateLayer(nEPLayer / 2);
322 const auto& senseLayer = geo.getSenseLayer(iSLayer);
323 const auto& fieldLayerIn = geo.getFieldLayer(iSLayer - 1);
324 const auto& fieldLayerOut = geo.getFieldLayer(iSLayer);
325
326 rmin_sensitive_left = epLayerBwd.getRmax();
327 rmax_sensitive_left = fieldLayerOut.getR();
328 zback_sensitive_left = senseLayer.getZbwd();
329 zfor_sensitive_left = epLayerBwd.getZfwd();
330
331 rmin_sensitive_middle = fieldLayerIn.getR();
332 rmax_sensitive_middle = fieldLayerOut.getR();
333 zback_sensitive_middle = epLayerBwd.getZfwd();
334 zfor_sensitive_middle = epLayerFwd.getZbwd();
335
336 rmin_sensitive_right = epLayerFwd.getRmax();
337 rmax_sensitive_right = fieldLayerOut.getR();
338 zback_sensitive_right = epLayerFwd.getZbwd();
339 zfor_sensitive_right = senseLayer.getZfwd();
340
341 } else if (iSLayer == 55) {
342
343 const auto& epLayerBwdIn = endplate.getEndPlateLayer(0);
344 const auto& epLayerBwdOut = endplate.getEndPlateLayer((nEPLayer / 2) - 1);
345 const auto& epLayerFwdIn = endplate.getEndPlateLayer(nEPLayer / 2);
346 const auto& epLayerFwdOut = endplate.getEndPlateLayer(nEPLayer - 1);
347 const auto& senseLayer = geo.getSenseLayer(iSLayer);
348 // const auto& fieldLayerIn = geo.getFieldLayer(iSLayer - 1);
349 int iSLayerMinus1 = iSLayer - 1; //avoid cpp-check warning
350 const auto& fieldLayerIn = geo.getFieldLayer(iSLayerMinus1); //avoid cpp-check warning
351 rmin_sensitive_left = epLayerBwdIn.getRmax();
352 rmax_sensitive_left = epLayerBwdOut.getRmax();
353 zback_sensitive_left = senseLayer.getZbwd();
354 zfor_sensitive_left = epLayerBwdIn.getZfwd();
355
356 rmin_sensitive_middle = fieldLayerIn.getR();
357 rmax_sensitive_middle = (geo.getOuterWall(0)).getRmin();
358 zback_sensitive_middle = epLayerBwdIn.getZfwd();
359 zfor_sensitive_middle = epLayerFwdIn.getZbwd();
360
361 rmin_sensitive_right = epLayerFwdIn.getRmax();
362 rmax_sensitive_right = epLayerFwdOut.getRmax();
363 zback_sensitive_right = epLayerFwdIn.getZbwd();
364 zfor_sensitive_right = senseLayer.getZfwd();
365
366 } else {
367 B2ERROR("Undefined sensitive layer : " << iSLayer);
368 continue;
369 }
370
371
372 // Check if build left sensitive tube
373 if ((zfor_sensitive_left - zback_sensitive_left) > length_feedthrough) {
374 // std::cout <<"left doif " << iSLayer <<" "<< zfor_sensitive_left - zback_sensitive_left << std::endl;
375 //==========================================================
376 // zback_sensitive_left
377 // |
378 // \|/
379 // _____________________
380 // | |// 1 // | |
381 // | |==ft====|2 | (ft = feedthrouth)
382 // |_______|____1___|__|
383 // |_______|___________|
384 // |
385 // \|/
386 // zfor_sensitive_left
387 //==========================================================
388
389 // Build a tube with metarial cdcMed for area 1
390 G4Tubs* leftTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_leftTube") % iSLayer).str().c_str(),
391 rmin_sensitive_left * CLHEP::cm,
392 rmax_sensitive_left * CLHEP::cm, length_feedthrough * CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
393 G4LogicalVolume* leftTube = new G4LogicalVolume(leftTubeShape, cdcMed,
394 (boost::format("logicalCDCLayer_%1%_leftTube") % iSLayer).str().c_str(), 0, 0, 0);
395 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zback_sensitive_left + length_feedthrough / 2.0)*CLHEP::cm), leftTube,
396 (boost::format("physicalCDCLayer_%1%_leftTube") % iSLayer).str().c_str(), m_logicalCDC, false, iSLayer);
397 // Build left sensitive tube (area 2)
398 G4Tubs* leftSensitiveTubeShape = new G4Tubs((boost::format("solidSD_CDCLayer_%1%_left") % iSLayer).str().c_str(),
399 rmin_sensitive_left * CLHEP::cm, rmax_sensitive_left * CLHEP::cm,
400 (zfor_sensitive_left - zback_sensitive_left - length_feedthrough)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
401 G4LogicalVolume* leftSensitiveTube = new G4LogicalVolume(leftSensitiveTubeShape, cdcMed,
402 (boost::format("logicalSD_CDCLayer_%1%_left") % iSLayer).str().c_str(), 0, 0, 0);
403 leftSensitiveTube->SetSensitiveDetector(m_sensitive);
404 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_left + zback_sensitive_left + length_feedthrough)*CLHEP::cm / 2.0),
405 leftSensitiveTube, (boost::format("physicalSD_CDCLayer_%1%_left") % iSLayer).str().c_str(), m_logicalCDC, false, iSLayer);
406 } else {
407 // std::cout <<"left doelse " << iSLayer << std::endl;
408 //==========================================================
409 // zback_sensitive_left
410 // |
411 // \|/
412 // _________________________
413 // | |//// 1 ////| 2 |
414 // | |======ft======== (ft = feedthrouth)
415 // |_______|____1______| 2 |
416 // |_______|___________|___|
417 // |
418 // \|/
419 // zfor_sensitive_left
420 //==========================================================
421
422 // Build a tube with metarial cdcMed for area 1
423 G4Tubs* leftTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_leftTube") % iSLayer).str().c_str(),
424 rmin_sensitive_left * CLHEP::cm,
425 rmax_sensitive_left * CLHEP::cm, (zfor_sensitive_left - zback_sensitive_left)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
426 G4LogicalVolume* leftTube = new G4LogicalVolume(leftTubeShape, cdcMed,
427 (boost::format("logicalCDCLayer_%1%_leftTube") % iSLayer).str().c_str(), 0, 0, 0);
428 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_left + zback_sensitive_left)*CLHEP::cm / 2.0), leftTube,
429 (boost::format("physicalCDCLayer_%1%_leftTube") % iSLayer).str().c_str(), m_logicalCDC, false, iSLayer);
430
431
432 // Build a tube with metarial cdcMed for area 2
433 G4Tubs* leftMidTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_leftMidTube") % iSLayer).str().c_str(),
434 rmin_sensitive_middle * CLHEP::cm, rmax_sensitive_middle * CLHEP::cm,
435 (length_feedthrough - zfor_sensitive_left + zback_sensitive_left)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
436 G4LogicalVolume* leftMidTube = new G4LogicalVolume(leftMidTubeShape, cdcMed,
437 (boost::format("logicalCDCLayer_%1%_leftMidTube") % iSLayer).str().c_str(), 0, 0, 0);
438
439 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (length_feedthrough + zfor_sensitive_left + zback_sensitive_left)*CLHEP::cm / 2.0),
440 leftMidTube, (boost::format("physicalCDCLayer_%1%_leftMidTube") % iSLayer).str().c_str(), m_logicalCDC, false, iSLayer);
441
442 // Reset zback_sensitive_middle
443 zback_sensitive_middle = length_feedthrough + zback_sensitive_left;
444 }
445
446 // Check if build right sensitive tube
447 if ((zfor_sensitive_right - zback_sensitive_right) > length_feedthrough) {
448 // std::cout <<"right doif" << iSLayer <<" "<< zfor_sensitive_right - zback_sensitive_right << std::endl;
449 //==========================================================
450 // zfor_sensitive_right
451 // |
452 // \|/
453 // _____________________________
454 // | | 1 |///////|
455 // | 2 |====ft=====|///////| (ft = feedthrouth)
456 // |_______|____1______|_______|
457 // |_______|___________|_______|
458 // |
459 // \|/
460 // zback_sensitive_right
461 //==========================================================
462
463 // Build a tube with metarial cdcMed for area 1
464 G4Tubs* rightTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_rightTube") % iSLayer).str().c_str(),
465 rmin_sensitive_right * CLHEP::cm, rmax_sensitive_right * CLHEP::cm, length_feedthrough * CLHEP::cm / 2.0, 0 * CLHEP::deg,
466 360.*CLHEP::deg);
467 G4LogicalVolume* rightTube = new G4LogicalVolume(rightTubeShape, cdcMed,
468 (boost::format("logicalCDCLayer_%1%_rightTube") % iSLayer).str().c_str(), 0, 0, 0);
469
470 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_right - length_feedthrough / 2.0)*CLHEP::cm), rightTube,
471 (boost::format("physicalCDCLayer_%1%_rightTube") % iSLayer).str().c_str(), m_logicalCDC, false, iSLayer);
472
473
474 // Build right sensitive tube (area 2)
475 G4Tubs* rightSensitiveTubeShape = new G4Tubs((boost::format("solidSD_CDCLayer_%1%_right") % iSLayer).str().c_str(),
476 rmin_sensitive_right * CLHEP::cm, rmax_sensitive_right * CLHEP::cm,
477 (zfor_sensitive_right - zback_sensitive_right - length_feedthrough)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
478 G4LogicalVolume* rightSensitiveTube = new G4LogicalVolume(rightSensitiveTubeShape, cdcMed,
479 (boost::format("logicalSD_CDCLayer_%1%_right") % iSLayer).str().c_str(), 0, 0, 0);
480 rightSensitiveTube->SetSensitiveDetector(m_sensitive);
481
482 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_right + zback_sensitive_right - length_feedthrough)*CLHEP::cm / 2.0),
483 rightSensitiveTube, (boost::format("physicalSD_CDCLayer_%1%_right") % iSLayer).str().c_str(), m_logicalCDC, false, iSLayer);
484
485 } else {
486 // std::cout <<"right doelse" << iSLayer << std::endl;
487 //==========================================================
488 // zfor_sensitive_right
489 // |
490 // \|/
491 // _____________________________
492 // | | 1 |///////|
493 // |============ft=====|///////| (ft = feedthrouth)
494 // | |____1______|_______|
495 // |_______|___________|_______|
496 // |
497 // \|/
498 // zback_sensitive_right
499 //==========================================================
500
501 // Build a tube with metarial cdcMed for area 1
502 G4Tubs* rightTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_rightTube") % iSLayer).str().c_str(),
503 rmin_sensitive_right * CLHEP::cm, rmax_sensitive_right * CLHEP::cm, (zfor_sensitive_right - zback_sensitive_right)*CLHEP::cm / 2.0,
504 0 * CLHEP::deg, 360.*CLHEP::deg);
505 G4LogicalVolume* rightTube = new G4LogicalVolume(rightTubeShape, cdcMed,
506 (boost::format("logicalCDCLayer_%1%_rightTube") % iSLayer).str().c_str(), 0, 0, 0);
507
508 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_right + zback_sensitive_right)*CLHEP::cm / 2.0), rightTube,
509 (boost::format("physicalCDCLayer_%1%_rightTube") % iSLayer).str().c_str(), m_logicalCDC, false, iSLayer);
510
511
512 // Build a tube with metarial cdcMed for area 2
513 G4Tubs* rightMidTubeShape = new G4Tubs((boost::format("solidCDCLayer_%1%_rightMidTube") % iSLayer).str().c_str(),
514 rmin_sensitive_middle * CLHEP::cm, rmax_sensitive_middle * CLHEP::cm,
515 (length_feedthrough - zfor_sensitive_right + zback_sensitive_right)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
516 G4LogicalVolume* rightMidTube = new G4LogicalVolume(rightMidTubeShape, cdcMed,
517 (boost::format("logicalCDCLayer_%1%_rightMidTube") % iSLayer).str().c_str(), 0, 0, 0);
518 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zback_sensitive_right - length_feedthrough + zfor_sensitive_right)*CLHEP::cm / 2.0),
519 rightMidTube, (boost::format("physicalCDCLayer_%1%_rightMidTube") % iSLayer).str().c_str(), m_logicalCDC, false, iSLayer);
520
521 // Reset zback_sensitive_middle
522 zfor_sensitive_middle = zfor_sensitive_right - length_feedthrough;
523 }
524
525
526 // Middle sensitive tube
527 G4Tubs* middleSensitiveTubeShape = new G4Tubs((boost::format("solidSD_CDCLayer_%1%_middle") % iSLayer).str().c_str(),
528 rmin_sensitive_middle * CLHEP::cm, rmax_sensitive_middle * CLHEP::cm,
529 (zfor_sensitive_middle - zback_sensitive_middle)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
530 G4LogicalVolume* middleSensitiveTube = new G4LogicalVolume(middleSensitiveTubeShape, cdcMedGas,
531 (boost::format("logicalSD_CDCLayer_%1%_middle") % iSLayer).str().c_str(), 0, 0, 0);
532 //hard-coded temporarily
533 //need to create an object per layer ??? to be checked later
534 G4UserLimits* uLimits = new G4UserLimits(8.5 * CLHEP::cm);
535 m_userLimits.push_back(uLimits);
536 middleSensitiveTube->SetUserLimits(uLimits);
537 middleSensitiveTube->SetSensitiveDetector(m_sensitive);
538
539 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfor_sensitive_middle + zback_sensitive_middle)*CLHEP::cm / 2.0), middleSensitiveTube,
540 (boost::format("physicalSD_CDCLayer_%1%_middle") % iSLayer).str().c_str(), m_logicalCDC, false, iSLayer);
541
542 // if (cdcgp.getMaterialDefinitionMode() == 2) {
543 if (gcp.getMaterialDefinitionMode() == 2) {
544 G4String sName = "sWire";
545 const G4int jc = 0;
546 B2Vector3D wb0 = cdcgp.wireBackwardPosition(iSLayer, jc);
547 // G4double rsense0 = wb0.Perp();
548 B2Vector3D wf0 = cdcgp.wireForwardPosition(iSLayer, jc);
549 G4double tAtZ0 = -wb0.Z() / (wf0.Z() - wb0.Z()); //t: param. along the wire
550 B2Vector3D wAtZ0 = wb0 + tAtZ0 * (wf0 - wb0);
551 //additional chop of wire; must be larger than 126um*(1/10), where 126um is the field wire diameter; 1/10: approx. stereo angle
552 const G4double epsl = 25.e-4; // (in cm);
553 G4double reductionBwd = (zback_sensitive_middle + epsl) / wb0.Z();
554 //chop the wire at zback_sensitive_middle for avoiding overlap; this is because the wire length defined by wb0 and wf0 is larger than the length of middle sensitive tube
555 wb0 = reductionBwd * (wb0 - wAtZ0) + wAtZ0;
556 //chop wire at zfor_sensitive_middle for avoiding overlap
557 G4double reductionFwd = (zfor_sensitive_middle - epsl) / wf0.Z();
558 wf0 = reductionFwd * (wf0 - wAtZ0) + wAtZ0;
559
560 const G4double wireHalfLength = 0.5 * (wf0 - wb0).Mag() * CLHEP::cm;
561 const G4double sWireRadius = 0.5 * cdcgp.senseWireDiameter() * CLHEP::cm;
562 // const G4double sWireRadius = 15.e-4 * CLHEP::cm;
563 G4Tubs* middleSensitiveSwireShape = new G4Tubs(sName, 0., sWireRadius, wireHalfLength, 0., 360. * CLHEP::deg);
564 G4LogicalVolume* middleSensitiveSwire = new G4LogicalVolume(middleSensitiveSwireShape, medTungsten, sName);
565 // middleSensitiveSwire->SetSensitiveDetector(m_sensitive);
566 middleSensitiveSwire->SetVisAttributes(m_VisAttributes.front()); // <- to speed up visualization
567
568 G4String fName = "fWire";
569 const G4double fWireRadius = 0.5 * cdcgp.fieldWireDiameter() * CLHEP::cm;
570 G4Tubs* middleSensitiveFwireShape = new G4Tubs(fName, 0., fWireRadius, wireHalfLength, 0., 360. * CLHEP::deg);
571 G4LogicalVolume* middleSensitiveFwire = new G4LogicalVolume(middleSensitiveFwireShape, medAluminum, fName);
572 // middleSensitiveFwire->SetSensitiveDetector(m_sensitive);
573 middleSensitiveFwire->SetVisAttributes(m_VisAttributes.front()); // <- to speed up visualization
574
575 const G4double diameter = cdcgp.fieldWireDiameter();
576
577 const G4int nCells = cdcgp.nWiresInLayer(iSLayer);
578 const G4double dphi = M_PI / nCells;
579 const B2Vector3D unitZ(0., 0., 1.);
580
581 for (int ic = 0; ic < nCells; ++ic) {
582 //define sense wire
583 B2Vector3D wb = cdcgp.wireBackwardPosition(iSLayer, ic);
584 B2Vector3D wf = cdcgp.wireForwardPosition(iSLayer, ic);
585 G4double tAtZ02 = -wb.Z() / (wf.Z() - wb.Z());
586 B2Vector3D wAtZ02 = wb + tAtZ02 * (wf - wb);
587 G4double reductionBwd2 = (zback_sensitive_middle + epsl) / wb.Z();
588 wb = reductionBwd2 * (wb - wAtZ02) + wAtZ02;
589 G4double reductionFwd2 = (zfor_sensitive_middle - epsl) / wf.Z();
590 wf = reductionFwd2 * (wf - wAtZ02) + wAtZ02;
591
592 G4double thetaYZ = -asin((wf - wb).Y() / (wf - wb).Mag());
593
594 B2Vector3D fMinusBInZX((wf - wb).X(), 0., (wf - wb).Z());
595 G4double thetaZX = asin((unitZ.Cross(fMinusBInZX)).Y() / fMinusBInZX.Mag());
596 G4RotationMatrix rotM;
597 // std::cout <<"deg,rad= " << CLHEP::deg <<" "<< CLHEP::rad << std::endl;
598 rotM.rotateX(thetaYZ * CLHEP::rad);
599 rotM.rotateY(thetaZX * CLHEP::rad);
600
601 G4ThreeVector xyz(0.5 * (wb.X() + wf.X()) * CLHEP::cm,
602 0.5 * (wb.Y() + wf.Y()) * CLHEP::cm, 0.);
603
604 // std::cout <<"0 x,y= " << xyz.getX() <<" "<< xyz.getY() << std::endl;
605 //Calling G4PVPlacement with G4Transform3D is convenient because it rotates the object instead of rotating the coordinate-frame; rotM is copied so it does not have to be created on heep by new.
606 new G4PVPlacement(G4Transform3D(rotM, xyz), middleSensitiveSwire, sName, middleSensitiveTube, false, ic);
607
608 //define field wire #1 (placed at the same phi but at the outer r boundary)
609 B2Vector3D wbF = wb;
610 G4double rF = rmax_sensitive_middle - 0.5 * diameter;
611 // std::cout <<"iSLayer,rF= " << iSLayer <<" "<< rF <<" "<< std::endl;
612 G4double phi = atan2(wbF.Y(), wbF.X());
613 wbF.SetX(rF * cos(phi));
614 wbF.SetY(rF * sin(phi));
615
616 B2Vector3D wfF = wf;
617 rF = rmax_sensitive_middle - 0.5 * diameter;
618 phi = atan2(wfF.Y(), wfF.X());
619 wfF.SetX(rF * cos(phi));
620 wfF.SetY(rF * sin(phi));
621
622 thetaYZ = -asin((wfF - wbF).Y() / (wfF - wbF).Mag());
623
624 fMinusBInZX = wfF - wbF;
625 fMinusBInZX.SetY(0.);
626 thetaZX = asin((unitZ.Cross(fMinusBInZX)).Y() / fMinusBInZX.Mag());
627
628 G4RotationMatrix rotM1;
629 rotM1.rotateX(thetaYZ * CLHEP::rad);
630 rotM1.rotateY(thetaZX * CLHEP::rad);
631
632 xyz.setX(0.5 * (wbF.X() + wfF.X()) * CLHEP::cm);
633 xyz.setY(0.5 * (wbF.Y() + wfF.Y()) * CLHEP::cm);
634
635 if (iSLayer != nSLayer - 1) {
636 // std::cout <<"1 x,y= " << xyz.getX() <<" "<< xyz.getY() << std::endl;
637 new G4PVPlacement(G4Transform3D(rotM1, xyz), middleSensitiveFwire, fName, middleSensitiveTube, false, ic);
638 }
639
640 //define field wire #2 (placed at the same radius but shifted by dphi)
641 wbF = wb;
642 rF = wbF.Perp();
643 phi = atan2(wbF.Y(), wbF.X());
644 wbF.SetX(rF * cos(phi + dphi));
645 wbF.SetY(rF * sin(phi + dphi));
646
647 wfF = wf;
648 rF = wfF.Perp();
649 phi = atan2(wfF.Y(), wfF.X());
650 wfF.SetX(rF * cos(phi + dphi));
651 wfF.SetY(rF * sin(phi + dphi));
652
653 thetaYZ = -asin((wfF - wbF).Y() / (wfF - wbF).Mag());
654
655 fMinusBInZX = wfF - wbF;
656 fMinusBInZX.SetY(0.);
657 thetaZX = asin((unitZ.Cross(fMinusBInZX)).Y() / fMinusBInZX.Mag());
658
659 G4RotationMatrix rotM2;
660 rotM2.rotateX(thetaYZ * CLHEP::rad);
661 rotM2.rotateY(thetaZX * CLHEP::rad);
662
663 xyz.setX(0.5 * (wbF.X() + wfF.X()) * CLHEP::cm);
664 xyz.setY(0.5 * (wbF.Y() + wfF.Y()) * CLHEP::cm);
665
666 // std::cout <<"2 x,y= " << xyz.getX() <<" "<< xyz.getY() << std::endl;
667 new G4PVPlacement(G4Transform3D(rotM2, xyz), middleSensitiveFwire, fName, middleSensitiveTube, false, ic + nCells);
668
669 //define field wire #3 (placed at the cell corner)
670 wbF = wb;
671 rF = rmax_sensitive_middle - 0.5 * diameter;
672 phi = atan2(wbF.Y(), wbF.X());
673 wbF.SetX(rF * cos(phi + dphi));
674 wbF.SetY(rF * sin(phi + dphi));
675
676 wfF = wf;
677 rF = rmax_sensitive_middle - 0.5 * diameter;
678 phi = atan2(wfF.Y(), wfF.X());
679 wfF.SetX(rF * cos(phi + dphi));
680 wfF.SetY(rF * sin(phi + dphi));
681
682 thetaYZ = -asin((wfF - wbF).Y() / (wfF - wbF).Mag());
683
684 fMinusBInZX = wfF - wbF;
685 fMinusBInZX.SetY(0.);
686 thetaZX = asin((unitZ.Cross(fMinusBInZX)).Y() / fMinusBInZX.Mag());
687
688 G4RotationMatrix rotM3;
689 rotM3.rotateX(thetaYZ * CLHEP::rad);
690 rotM3.rotateY(thetaZX * CLHEP::rad);
691
692 xyz.setX(0.5 * (wbF.X() + wfF.X()) * CLHEP::cm);
693 xyz.setY(0.5 * (wbF.Y() + wfF.Y()) * CLHEP::cm);
694
695 if (iSLayer != nSLayer - 1) {
696 new G4PVPlacement(G4Transform3D(rotM3, xyz), middleSensitiveFwire, fName, middleSensitiveTube, false, ic + 2 * nCells);
697 }
698 } // end of wire loop
699 } // end of wire definitions
700
701 }
702 //
703 // Endplates.
704 //
705
706 m_VisAttributes.push_back(new G4VisAttributes(true, G4Colour(1., 1., 0.)));
707 for (const auto& endplate : geo.getEndPlates()) {
708 for (const auto& epLayer : endplate.getEndPlateLayers()) {
709 const int iEPLayer = epLayer.getILayer();
710 const string name = epLayer.getName();
711 const double rmin = epLayer.getRmin();
712 const double rmax = epLayer.getRmax();
713 const double zbwd = epLayer.getZbwd();
714 const double zfwd = epLayer.getZfwd();
715 const double length = (zfwd - zbwd) / 2.0;
716
717 G4Tubs* tube = new G4Tubs("solidCDCEndplate" + name, rmin * CLHEP::cm,
718 rmax * CLHEP::cm, length * CLHEP::cm, 0 * CLHEP::deg, 360.*CLHEP::deg);
719 G4LogicalVolume* logical = new G4LogicalVolume(tube, Materials::get("G4_Al"),
720 "logicalCDCEndplate" + name, 0, 0);
721 logical->SetVisAttributes(m_VisAttributes.back());
722 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (zfwd + zbwd)*CLHEP::cm / 2.0), logical,
723 "physicalCDCEndplate" + name, m_logicalCDC, false, iEPLayer);
724
725 }
726 }
727
728
729 // Construct electronics boards
730 for (const auto& frontend : geo.getFrontends()) {
731
732 const int iEB = frontend.getId();
733 const double ebInnerR = frontend.getRmin();
734 const double ebOuterR = frontend.getRmax();
735 const double ebBZ = frontend.getZbwd();
736 const double ebFZ = frontend.getZfwd();
737
738 G4Tubs* ebTubeShape = new G4Tubs((boost::format("solidSD_ElectronicsBoard_Layer%1%") % iEB).str().c_str(), ebInnerR * CLHEP::cm,
739 ebOuterR * CLHEP::cm, (ebFZ - ebBZ)*CLHEP::cm / 2.0, 0 * CLHEP::deg, 360.*CLHEP::deg);
740
741 G4LogicalVolume* ebTube = new G4LogicalVolume(ebTubeShape, medNEMA_G10_Plate,
742 (boost::format("logicalSD_ElectronicsBoard_Layer%1%") % iEB).str().c_str(), 0, 0, 0);
743 if (!m_bkgsensitive) m_bkgsensitive = new BkgSensitiveDetector("CDC", iEB);
744 ebTube->SetSensitiveDetector(m_bkgsensitive);
745 ebTube->SetVisAttributes(m_VisAttributes.back());
746 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (ebFZ + ebBZ)*CLHEP::cm / 2.0), ebTube,
747 (boost::format("physicalSD_ElectronicsBoard_Layer%1%") % iEB).str().c_str(), m_logicalCDC, false, iEB);
748 }
749
750 //
751 // Construct neutron shield.
752 //
754
755 //
756 // construct covers.
757 //
758 createCovers(geo);
759
760 //
761 // construct covers.
762 //
763 createCover2s(geo);
764
765 //
766 // Construct ribs.
767 //
768 for (const auto& rib : geo.getRibs()) {
769
770 const int id = rib.getId();
771 const double length = rib.getLength();
772 const double width = rib.getWidth();
773 const double thick = rib.getThick();
774 const double rotx = rib.getRotX();
775 const double roty = rib.getRotY();
776 const double rotz = rib.getRotZ();
777 const double x = rib.getX();
778 const double y = rib.getY();
779 const double z = rib.getZ();
780 const int offset = rib.getOffset();
781 const int ndiv = rib.getNDiv();
782
783 const string solidName = "solidRib" + to_string(id);
784 const string logicalName = "logicalRib" + to_string(id);
785 G4Box* boxShape = new G4Box(solidName, 0.5 * length * CLHEP::cm,
786 0.5 * width * CLHEP::cm,
787 0.5 * thick * CLHEP::cm);
788
789 const double rmax = 0.5 * length;
790 const double rmin = max((rmax - thick), 0.);
791 G4Tubs* tubeShape = new G4Tubs(solidName,
792 rmin * CLHEP::cm,
793 rmax * CLHEP::cm,
794 0.5 * width * CLHEP::cm,
795 0.,
796 360. * CLHEP::deg);
797
798 //G4LogicalVolume* logicalV = new G4LogicalVolume(boxShape, medAluminum,
799 // logicalName, 0, 0, 0);
800 // ID dependent material definition, Aluminum is default
801 G4LogicalVolume* logicalV = new G4LogicalVolume(boxShape, medAluminum, logicalName, 0, 0, 0);
802 if (id > 39 && id < 78) // Cu
803 logicalV = new G4LogicalVolume(boxShape, medCopper, logicalName, 0, 0, 0);
804 if ((id > 77 && id < 94) || (id > 131 && id < 146)) // G10
805 logicalV = new G4LogicalVolume(boxShape, medNEMA_G10_Plate, logicalName, 0, 0, 0);
806 if (id > 93 && id < 110) // Cu
807 logicalV = new G4LogicalVolume(tubeShape, medCopper, logicalName, 0, 0, 0);
808 if (id > 109 && id < 126) // H2O
809 logicalV = new G4LogicalVolume(tubeShape, medH2O, logicalName, 0, 0, 0);
810 if (id > 127 && id < 132) // HV
811 logicalV = new G4LogicalVolume(boxShape, medHV, logicalName, 0, 0, 0);
812 /*if( id > 145 && id < 149 )// Fiber
813 logicalV = new G4LogicalVolume(boxShape, medFiber, logicalName, 0, 0, 0);
814 if( id > 148 && id < 158 )// Fiber
815 logicalV = new G4LogicalVolume(boxShape, medCAT7, logicalName, 0, 0, 0);
816 if( id > 157 && id < 164 )// Fiber
817 logicalV = new G4LogicalVolume(boxShape, medTRG, logicalName, 0, 0, 0);*/
818
819 logicalV->SetVisAttributes(m_VisAttributes.back());
820
821 const double phi = 360.0 / ndiv;
822
823 G4RotationMatrix rot = G4RotationMatrix();
824 double dz = thick;
825 if (id > 93 && id < 126) dz = 0;
826
827 G4ThreeVector arm(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - dz * CLHEP::cm / 2.0);
828 rot.rotateX(rotx);
829 rot.rotateY(roty);
830 rot.rotateZ(rotz);
831 if (offset) {
832 rot.rotateZ(0.5 * phi * CLHEP::deg);
833 arm.rotateZ(0.5 * phi * CLHEP::deg);
834 }
835 for (int i = 0; i < ndiv; ++i) {
836 const string physicalName = "physicalRib_" + to_string(id) + " " + to_string(i);
837 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
838 physicalName.c_str(), m_logicalCDC, false, id);
839 rot.rotateZ(phi * CLHEP::deg);
840 arm.rotateZ(phi * CLHEP::deg);
841 }
842
843 }
844
845 //
846 // Construct rib2s.
847 //
848 for (const auto& rib2 : geo.getRib2s()) {
849
850 const int id = rib2.getId();
851 const double length = rib2.getLength();
852 const double width = rib2.getWidth();
853 const double thick = rib2.getThick();
854 const double width2 = rib2.getWidth2();
855 const double thick2 = rib2.getThick2();
856 const double rotx = rib2.getRotX();
857 const double roty = rib2.getRotY();
858 const double rotz = rib2.getRotZ();
859 const double x = rib2.getX();
860 const double y = rib2.getY();
861 const double z = rib2.getZ();
862 const int ndiv = rib2.getNDiv();
863
864 const string solidName = "solidRib2" + to_string(id);
865 const string logicalName = "logicalRib2" + to_string(id);
866 G4Trd* trdShape = new G4Trd(solidName,
867 0.5 * thick * CLHEP::cm,
868 0.5 * thick2 * CLHEP::cm,
869 0.5 * width * CLHEP::cm,
870 0.5 * width2 * CLHEP::cm,
871 0.5 * length * CLHEP::cm);
872
873 G4LogicalVolume* logicalV = new G4LogicalVolume(trdShape, medAluminum, logicalName, 0, 0, 0);
874
875 if (id > 0)
876 logicalV = new G4LogicalVolume(trdShape, medCopper, logicalName, 0, 0, 0);
877
878 logicalV->SetVisAttributes(m_VisAttributes.back());
879
880 const double phi = 360.0 / ndiv;
881
882 G4RotationMatrix rot = G4RotationMatrix();
883 G4ThreeVector arm(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - thick * CLHEP::cm / 2.0);
884
885 rot.rotateX(rotx);
886 rot.rotateY(roty);
887 rot.rotateZ(rotz);
888 for (int i = 0; i < ndiv; ++i) {
889 const string physicalName = "physicalRib2_" + to_string(id) + " " + to_string(i);
890 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
891 physicalName.c_str(), m_logicalCDC, false, id);
892 rot.rotateZ(phi * CLHEP::deg);
893 arm.rotateZ(phi * CLHEP::deg);
894 }
895
896 }
897
898 //
899 // Construct rib3s.
900 //
901 for (const auto& rib3 : geo.getRib3s()) {
902
903 const int id = rib3.getId();
904 const double length = rib3.getLength();
905 const double width = rib3.getWidth();
906 const double thick = rib3.getThick();
907 const double r = rib3.getR();
908 const double x = rib3.getX();
909 const double y = rib3.getY();
910 const double z = rib3.getZ();
911 const double rx = rib3.getRx();
912 const double ry = rib3.getRy();
913 const double rz = rib3.getRz();
914 const int offset = rib3.getOffset();
915 const int ndiv = rib3.getNDiv();
916
917 const string logicalName = "logicalRib3" + to_string(id);
918 G4VSolid* boxShape = new G4Box("Block",
919 0.5 * length * CLHEP::cm,
920 0.5 * width * CLHEP::cm,
921 0.5 * thick * CLHEP::cm);
922 G4VSolid* tubeShape = new G4Tubs("Hole",
923 0.,
924 r * CLHEP::cm,
925 width * CLHEP::cm,
926 0. * CLHEP::deg,
927 360. * CLHEP::deg);
928
929 G4RotationMatrix rotsub = G4RotationMatrix();
930 rotsub.rotateX(90. * CLHEP::deg);
931 G4ThreeVector trnsub(rx * CLHEP::cm - x * CLHEP::cm, ry * CLHEP::cm - y * CLHEP::cm,
932 rz * CLHEP::cm - z * CLHEP::cm + 0.5 * thick * CLHEP::cm);
933 G4VSolid* coolingBlock = new G4SubtractionSolid("Block-Hole",
934 boxShape,
935 tubeShape,
936 G4Transform3D(rotsub,
937 trnsub));
938
939 G4LogicalVolume* logicalV = new G4LogicalVolume(coolingBlock, medCopper, logicalName, 0, 0, 0);
940
941 logicalV->SetVisAttributes(m_VisAttributes.back());
942
943 const double phi = 360.0 / ndiv;
944
945 G4RotationMatrix rot = G4RotationMatrix();
946 G4ThreeVector arm(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - thick * CLHEP::cm / 2.0);
947
948 if (offset) {
949 rot.rotateZ(0.5 * phi * CLHEP::deg);
950 arm.rotateZ(0.5 * phi * CLHEP::deg);
951 }
952 for (int i = 0; i < ndiv; ++i) {
953 const string physicalName = "physicalRib3_" + to_string(id) + " " + to_string(i);
954 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
955 physicalName.c_str(), m_logicalCDC, false, id);
956 rot.rotateZ(phi * CLHEP::deg);
957 arm.rotateZ(phi * CLHEP::deg);
958 }
959
960 }
961
962 //
963 // Construct rib4s.
964 //
965 for (const auto& rib4 : geo.getRib4s()) {
966
967 const int id = rib4.getId();
968 const double length = rib4.getLength();
969 const double width = rib4.getWidth();
970 const double thick = rib4.getThick();
971 const double length2 = rib4.getLength2();
972 const double width2 = rib4.getWidth2();
973 const double thick2 = rib4.getThick2();
974 const double x = rib4.getX();
975 const double y = rib4.getY();
976 const double z = rib4.getZ();
977 const double x2 = rib4.getX2();
978 const double y2 = rib4.getY2();
979 const double z2 = rib4.getZ2();
980 const int offset = rib4.getOffset();
981 const int ndiv = rib4.getNDiv();
982
983 const string logicalName = "logicalRib4" + to_string(id);
984 G4VSolid* baseShape = new G4Box("Base",
985 0.5 * length * CLHEP::cm,
986 0.5 * width * CLHEP::cm,
987 0.5 * thick * CLHEP::cm);
988 G4VSolid* sqShape = new G4Box("Sq",
989 0.5 * length2 * CLHEP::cm,
990 0.5 * width2 * CLHEP::cm,
991 0.5 * thick2 * CLHEP::cm);
992
993 G4RotationMatrix rotsub = G4RotationMatrix();
994 double dzc = (z2 - thick2 / 2.) - (z - thick / 2.);
995 G4ThreeVector trnsub(x2 * CLHEP::cm - x * CLHEP::cm,
996 y2 * CLHEP::cm - y * CLHEP::cm,
997 dzc * CLHEP::cm);
998 G4VSolid* sqHoleBase = new G4SubtractionSolid("Box-Sq",
999 baseShape,
1000 sqShape,
1001 G4Transform3D(rotsub,
1002 trnsub)
1003 );
1004
1005 G4LogicalVolume* logicalV = new G4LogicalVolume(sqHoleBase, medCopper, logicalName, 0, 0, 0);
1006 if (id < 19) {
1007 logicalV = new G4LogicalVolume(sqHoleBase, medNEMA_G10_Plate, logicalName, 0, 0, 0);
1008 BkgSensitiveDetector* sensitiveDetector =
1009 new BkgSensitiveDetector("CDC", 2000 + id);
1010 logicalV->SetSensitiveDetector(sensitiveDetector);
1011 m_BkgSensitiveRib4.push_back(sensitiveDetector);
1012 }
1013
1014 logicalV->SetVisAttributes(m_VisAttributes.back());
1015
1016 const double phi = 360.0 / ndiv;
1017
1018 G4RotationMatrix rot = G4RotationMatrix();
1019 G4ThreeVector arm(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - thick * CLHEP::cm / 2.0);
1020
1021 if (offset) {
1022 rot.rotateZ(0.5 * phi * CLHEP::deg);
1023 arm.rotateZ(0.5 * phi * CLHEP::deg);
1024 }
1025 for (int i = 0; i < ndiv; ++i) {
1026 const string physicalName = "physicalRib4_" + to_string(id) + " " + to_string(i);
1027 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
1028 physicalName.c_str(), m_logicalCDC, false, id);
1029 rot.rotateZ(phi * CLHEP::deg);
1030 arm.rotateZ(phi * CLHEP::deg);
1031 }
1032
1033 }
1034 //
1035 // Construct rib5s.
1036 //
1037 for (const auto& rib5 : geo.getRib5s()) {
1038
1039 const int id = rib5.getId();
1040 const double dr = rib5.getDr();
1041 const double dz = rib5.getDz();
1042 const double width = rib5.getWidth();
1043 const double thick = rib5.getThick();
1044 const double rin = rib5.getRin();
1045 const double x = rib5.getX();
1046 const double y = rib5.getY();
1047 const double z = rib5.getZ();
1048 const double rotx = rib5.getRotx();
1049 const double roty = rib5.getRoty();
1050 const double rotz = rib5.getRotz();
1051 const int offset = rib5.getOffset();
1052 const int ndiv = rib5.getNDiv();
1053
1054 const string solidName = "solidRib5" + to_string(id);
1055 const string logicalName = "logicalRib5" + to_string(id);
1056
1057 const double rmax = rin + thick;
1058 const double rmin = rin;
1059 const double dphi = 2. * atan2(dz, dr);
1060 const double ddphi = thick * tan(dphi) / rin;
1061 const double ddphi2 = width / 2. * width / 2. / (x + dr) / rin;
1062 const double cphi = dphi - ddphi - ddphi2;
1063 G4Tubs* tubeShape = new G4Tubs(solidName,
1064 rmin * CLHEP::cm,
1065 rmax * CLHEP::cm,
1066 0.5 * width * CLHEP::cm,
1067 0.,
1068 cphi);
1069
1070 G4LogicalVolume* logicalV = new G4LogicalVolume(tubeShape, medAluminum, logicalName, 0, 0, 0);
1071
1072 logicalV->SetVisAttributes(m_VisAttributes.back());
1073
1074 const double phi = 360.0 / ndiv;
1075
1076 G4RotationMatrix rot = G4RotationMatrix();
1077
1078 //G4ThreeVector arm(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - thick * CLHEP::cm / 2.0);
1079 G4ThreeVector arm(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - rin * CLHEP::cm - thick * CLHEP::cm);
1080 rot.rotateX(rotx);
1081 rot.rotateY(roty);
1082 rot.rotateZ(rotz);
1083 if (offset) {
1084 rot.rotateZ(0.5 * phi * CLHEP::deg);
1085 arm.rotateZ(0.5 * phi * CLHEP::deg);
1086 }
1087 for (int i = 0; i < ndiv; ++i) {
1088 const string physicalName = "physicalRib5_" + to_string(id) + " " + to_string(i);
1089 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
1090 physicalName.c_str(), m_logicalCDC, false, id);
1091 rot.rotateZ(phi * CLHEP::deg);
1092 arm.rotateZ(phi * CLHEP::deg);
1093 }
1094
1095 }
1096
1097 //Create B-field mapper (here tentatively)
1098 createMapper(topVolume);
1099 }
1100
1101
1103 {
1104
1105 G4Material* C2H4 = geometry::Materials::get("G4_POLYETHYLENE");
1106 G4Material* elB = geometry::Materials::get("G4_B");
1107
1108 // 5% borated polyethylene = SWX201
1109 // http://www.deqtech.com/Shieldwerx/Products/swx201hd.htm
1110 G4Material* boratedpoly05 = new G4Material("BoratedPoly05", 1.06 * CLHEP::g / CLHEP::cm3, 2);
1111 boratedpoly05->AddMaterial(elB, 0.05);
1112 boratedpoly05->AddMaterial(C2H4, 0.95);
1113 // 30% borated polyethylene = SWX210
1114 G4Material* boratedpoly30 = new G4Material("BoratedPoly30", 1.19 * CLHEP::g / CLHEP::cm3, 2);
1115 boratedpoly30->AddMaterial(elB, 0.30);
1116 boratedpoly30->AddMaterial(C2H4, 0.70);
1117
1118 G4Material* shieldMat = C2H4;
1119
1120 const int nShields = content.getNumberNodes("Shields/Shield");
1121
1122 for (int iShield = 0; iShield < nShields; ++iShield) {
1123 GearDir shieldContent(content);
1124 shieldContent.append((boost::format("/Shields/Shield[%1%]/") % (iShield + 1)).str());
1125 const string sShieldID = shieldContent.getString("@id");
1126 const int shieldID = atoi(sShieldID.c_str());
1127 // const string shieldName = shieldContent.getString("Name");
1128 const double shieldInnerR1 = shieldContent.getLength("InnerR1");
1129 const double shieldInnerR2 = shieldContent.getLength("InnerR2");
1130 const double shieldOuterR1 = shieldContent.getLength("OuterR1");
1131 const double shieldOuterR2 = shieldContent.getLength("OuterR2");
1132 const double shieldThick = shieldContent.getLength("Thickness");
1133 const double shieldPosZ = shieldContent.getLength("PosZ");
1134
1135 G4Cons* shieldConsShape = new G4Cons((boost::format("solidShield%1%") % shieldID).str().c_str(),
1136 shieldInnerR1 * CLHEP::cm, shieldOuterR1 * CLHEP::cm,
1137 shieldInnerR2 * CLHEP::cm, shieldOuterR2 * CLHEP::cm,
1138 shieldThick * CLHEP::cm / 2.0,
1139 0.*CLHEP::deg, 360.*CLHEP::deg);
1140
1141 G4LogicalVolume* shieldCons = new G4LogicalVolume(shieldConsShape, shieldMat,
1142 (boost::format("logicalShield%1%") % shieldID).str().c_str(),
1143 0, 0, 0);
1144 shieldCons->SetVisAttributes(m_VisAttributes.back());
1145 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (shieldPosZ - shieldThick / 2.0) * CLHEP::cm), shieldCons,
1146 (boost::format("physicalShield%1%") % shieldID).str().c_str(), m_logicalCDC, false, 0);
1147
1148 }
1149
1150 }
1151
1152
1154 {
1155
1156 G4Material* C2H4 = geometry::Materials::get("G4_POLYETHYLENE");
1157 G4Material* shieldMat = C2H4;
1158
1159 for (const auto& shield : geom.getNeutronShields()) {
1160 const int shieldID = shield.getId();
1161 const double shieldInnerR1 = shield.getRmin1();
1162 const double shieldInnerR2 = shield.getRmin2();
1163 const double shieldOuterR1 = shield.getRmax1();
1164 const double shieldOuterR2 = shield.getRmax2();
1165 const double shieldThick = shield.getThick();
1166 const double shieldPosZ = shield.getZ();
1167
1168 G4Cons* shieldConsShape = new G4Cons("solidShield" + to_string(shieldID),
1169 shieldInnerR1 * CLHEP::cm, shieldOuterR1 * CLHEP::cm,
1170 shieldInnerR2 * CLHEP::cm, shieldOuterR2 * CLHEP::cm,
1171 shieldThick * CLHEP::cm / 2.0,
1172 0.*CLHEP::deg, 360.*CLHEP::deg);
1173
1174 G4LogicalVolume* shieldCons = new G4LogicalVolume(shieldConsShape, shieldMat, "logicalShield" + to_string(shieldID),
1175 0, 0, 0);
1176 shieldCons->SetVisAttributes(m_VisAttributes.back());
1177 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, (shieldPosZ - shieldThick / 2.0) * CLHEP::cm), shieldCons,
1178 "physicalShield" + to_string(shieldID), m_logicalCDC, false, 0);
1179
1180 }
1181
1182 }
1183
1185 {
1186 string Aluminum = content.getString("Aluminum");
1187 G4Material* medAluminum = geometry::Materials::get(Aluminum);
1188 G4Material* medNEMA_G10_Plate = geometry::Materials::get("NEMA_G10_Plate");
1189 G4double density = 1.000 * CLHEP::g / CLHEP::cm3;
1190 G4double a = 1.01 * CLHEP::g / CLHEP::mole;
1191 G4Element* elH = new G4Element("Hydrogen", "H", 1., a);
1192 a = 16.00 * CLHEP::g / CLHEP::mole;
1193 G4Element* elO = new G4Element("Oxygen", "O", 8., a);
1194 G4Material* medH2O = new G4Material("Water", density, 2);
1195 medH2O->AddElement(elH, 2);
1196 medH2O->AddElement(elO, 1);
1197 G4Material* medCopper = geometry::Materials::get("Cu");
1198 G4Material* medLV = geometry::Materials::get("CDCLVCable");
1199 G4Material* medFiber = geometry::Materials::get("CDCOpticalFiber");
1200 G4Material* medCAT7 = geometry::Materials::get("CDCCAT7");
1201 G4Material* medTRG = geometry::Materials::get("CDCOpticalFiberTRG");
1202 G4Material* medHV = geometry::Materials::get("CDCHVCable");
1203
1204 m_VisAttributes.push_back(new G4VisAttributes(true, G4Colour(0., 1., 0.)));
1205 const int nCover = content.getNumberNodes("Covers/Cover");
1206 for (int iCover = 0; iCover < nCover; ++iCover) {
1207 GearDir coverContent(content);
1208 coverContent.append((boost::format("/Covers/Cover[%1%]/") % (iCover + 1)).str());
1209 const string scoverID = coverContent.getString("@id");
1210 const int coverID = atoi(scoverID.c_str());
1211 const string coverName = coverContent.getString("Name");
1212 const double coverInnerR1 = coverContent.getLength("InnerR1");
1213 const double coverInnerR2 = coverContent.getLength("InnerR2");
1214 const double coverOuterR1 = coverContent.getLength("OuterR1");
1215 const double coverOuterR2 = coverContent.getLength("OuterR2");
1216 const double coverThick = coverContent.getLength("Thickness");
1217 const double coverPosZ = coverContent.getLength("PosZ");
1218
1219 const double rmin1 = coverInnerR1;
1220 const double rmax1 = coverOuterR1;
1221 const double rmin2 = coverInnerR2;
1222 const double rmax2 = coverOuterR2;
1223
1224 /*
1225 if (coverID == 7 || coverID == 10) {
1226 createCone(rmin1, rmax1, rmin2, rmax2, coverThick, coverPosZ, coverID, medAluminum, coverName);
1227 } else {
1228 createTube(rmin1, rmax1, coverThick, coverPosZ, coverID, medAluminum, coverName);
1229
1230 }*/
1231 // ID dependent material definition
1232 if (coverID < 23) {
1233 if (coverID == 7 || coverID == 10) {// cones
1234 createCone(rmin1, rmax1, rmin2, rmax2, coverThick, coverPosZ, coverID, medAluminum, coverName);
1235 } else {// covers
1236 createTube(rmin1, rmax1, coverThick, coverPosZ, coverID, medAluminum, coverName);
1237 }
1238 }
1239 if (coverID > 22 && coverID < 29)// cooling plate
1240 createTube(rmin1, rmax1, coverThick, coverPosZ, coverID, medCopper, coverName);
1241 if (coverID > 28 && coverID < 35)// cooling Pipe
1242 createTorus(rmin1, rmax1, coverThick, coverPosZ, coverID, medCopper, coverName);
1243 if (coverID > 34 && coverID < 41)// cooling water
1244 createTorus(rmin1, rmax1, coverThick, coverPosZ, coverID, medH2O, coverName);
1245 if (coverID == 45 || coverID == 46)
1246 createTube(rmin1, rmax1, coverThick, coverPosZ, coverID, medLV, coverName);
1247 if (coverID == 47 || coverID == 48)
1248 createTube(rmin1, rmax1, coverThick, coverPosZ, coverID, medFiber, coverName);
1249 if (coverID == 49 || coverID == 50)
1250 createTube(rmin1, rmax1, coverThick, coverPosZ, coverID, medCAT7, coverName);
1251 if (coverID == 51 || coverID == 52)
1252 createTube(rmin1, rmax1, coverThick, coverPosZ, coverID, medTRG, coverName);
1253 if (coverID == 53)
1254 createTube(rmin1, rmax1, coverThick, coverPosZ, coverID, medHV, coverName);
1255 }
1256
1257 const int nCover2 = content.getNumberNodes("Covers/Cover2");
1258 for (int iCover2 = 0; iCover2 < nCover2; ++iCover2) {
1259 GearDir cover2Content(content);
1260 cover2Content.append((boost::format("/Cover2s/Cover2[%1%]/") % (iCover2 + 1)).str());
1261 const string scover2ID = cover2Content.getString("@id");
1262 const int cover2ID = atoi(scover2ID.c_str());
1263 const string cover2Name = cover2Content.getString("Name");
1264 const double cover2InnerR = cover2Content.getLength("InnerR");
1265 const double cover2OuterR = cover2Content.getLength("OuterR");
1266 const double cover2StartPhi = cover2Content.getLength("StartPhi");
1267 const double cover2DeltaPhi = cover2Content.getLength("DeltaPhi");
1268 const double cover2Thick = cover2Content.getLength("Thickness");
1269 const double cover2PosZ = cover2Content.getLength("PosZ");
1270
1271 if (cover2ID < 11)
1272 createTube2(cover2InnerR, cover2OuterR, cover2StartPhi, cover2DeltaPhi, cover2Thick, cover2PosZ, cover2ID, medHV, cover2Name);
1273 if (cover2ID > 10 && cover2ID < 14)
1274 createTube2(cover2InnerR, cover2OuterR, cover2StartPhi, cover2DeltaPhi, cover2Thick, cover2PosZ, cover2ID, medFiber, cover2Name);
1275 if (cover2ID > 13 && cover2ID < 23)
1276 createTube2(cover2InnerR, cover2OuterR, cover2StartPhi, cover2DeltaPhi, cover2Thick, cover2PosZ, cover2ID, medCAT7, cover2Name);
1277 if (cover2ID > 22 && cover2ID < 29)
1278 createTube2(cover2InnerR, cover2OuterR, cover2StartPhi, cover2DeltaPhi, cover2Thick, cover2PosZ, cover2ID, medTRG, cover2Name);
1279 }
1280
1281 const int nRibs = content.getNumberNodes("Covers/Rib");
1282 for (int iRib = 0; iRib < nRibs; ++iRib) {
1283 GearDir ribContent(content);
1284 ribContent.append((boost::format("/Covers/Rib[%1%]/") % (iRib + 1)).str());
1285 const string sribID = ribContent.getString("@id");
1286 const int ribID = atoi(sribID.c_str());
1287 // const string ribName = ribContent.getString("Name");
1288 const double length = ribContent.getLength("Length");
1289 const double width = ribContent.getLength("Width");
1290 const double thick = ribContent.getLength("Thickness");
1291 const double rotX = ribContent.getLength("RotX");
1292 const double rotY = ribContent.getLength("RotY");
1293 const double rotZ = ribContent.getLength("RotZ");
1294 const double cX = ribContent.getLength("PosX");
1295 const double cY = ribContent.getLength("PosY");
1296 const double cZ = ribContent.getLength("PosZ");
1297 const int offset = atoi((ribContent.getString("Offset")).c_str());
1298 const int number = atoi((ribContent.getString("NDiv")).c_str());
1299
1300 const string solidName = "solidRib" + to_string(ribID);
1301 const string logicalName = "logicalRib" + to_string(ribID);
1302 G4Box* boxShape = new G4Box(solidName, 0.5 * length * CLHEP::cm,
1303 0.5 * width * CLHEP::cm,
1304 0.5 * thick * CLHEP::cm);
1305 const double rmax = 0.5 * length;
1306 const double rmin = max((rmax - thick), 0.);
1307 G4Tubs* tubeShape = new G4Tubs(solidName,
1308 rmin * CLHEP::cm,
1309 rmax * CLHEP::cm,
1310 0.5 * width * CLHEP::cm,
1311 0.,
1312 360. * CLHEP::deg);
1313
1314 //G4LogicalVolume* logicalV = new G4LogicalVolume(boxShape, medAluminum,
1315 // logicalName, 0, 0, 0);
1316 // ID dependent material definition Aluminum is default
1317 G4LogicalVolume* logicalV = new G4LogicalVolume(boxShape, medAluminum, logicalName, 0, 0, 0);
1318 if (ribID > 39 && ribID < 78) // Cu box
1319 logicalV = new G4LogicalVolume(boxShape, medCopper, logicalName, 0, 0, 0);
1320 if ((ribID > 77 && ribID < 94) || (ribID > 131 && ribID < 146)) // G10 box
1321 logicalV = new G4LogicalVolume(boxShape, medNEMA_G10_Plate, logicalName, 0, 0, 0);
1322 if (ribID > 93 && ribID < 110) // Cu tube
1323 logicalV = new G4LogicalVolume(tubeShape, medCopper, logicalName, 0, 0, 0);
1324 if (ribID > 109 && ribID < 126) // H2O tube (rmin = 0)
1325 logicalV = new G4LogicalVolume(tubeShape, medH2O, logicalName, 0, 0, 0);
1326 if (ribID > 127 && ribID < 132) // HV bundle
1327 logicalV = new G4LogicalVolume(boxShape, medHV, logicalName, 0, 0, 0);
1328 /*if( ribID > 145 && ribID < 149 )// Fiber box
1329 logicalV = new G4LogicalVolume(boxShape, medFiber, logicalName, 0, 0, 0);
1330 if( ribID > 148 && ribID < 158 )// Fiber box
1331 logicalV = new G4LogicalVolume(boxShape, medCAT7, logicalName, 0, 0, 0);
1332 if( ribID > 157 && ribID < 164 )// Fiber box
1333 logicalV = new G4LogicalVolume(boxShape, medTRG, logicalName, 0, 0, 0);*/
1334
1335 logicalV->SetVisAttributes(m_VisAttributes.back());
1336
1337 const double phi = 360.0 / number;
1338
1339 G4RotationMatrix rot = G4RotationMatrix();
1340
1341 double dz = thick;
1342 if (ribID > 93 && ribID < 126) dz = 0;
1343 G4ThreeVector arm(cX * CLHEP::cm, cY * CLHEP::cm, cZ * CLHEP::cm - dz * CLHEP::cm / 2.0);
1344
1345 rot.rotateX(rotX);
1346 rot.rotateY(rotY);
1347 rot.rotateZ(rotZ);
1348 if (offset) {
1349 rot.rotateZ(0.5 * phi * CLHEP::deg);
1350 arm.rotateZ(0.5 * phi * CLHEP::deg);
1351 }
1352 for (int i = 0; i < number; ++i) {
1353 const string physicalName = "physicalRib_" + to_string(ribID) + " " + to_string(i);
1354 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
1355 physicalName.c_str(), m_logicalCDC, false, ribID);
1356 rot.rotateZ(phi * CLHEP::deg);
1357 arm.rotateZ(phi * CLHEP::deg);
1358 }
1359
1360 }// rib
1361
1362 const int nRib2s = content.getNumberNodes("Covers/Rib2");
1363 for (int iRib2 = 0; iRib2 < nRib2s; ++iRib2) {
1364 GearDir rib2Content(content);
1365 rib2Content.append((boost::format("/Covers/Rib2[%1%]/") % (iRib2 + 1)).str());
1366 const string srib2ID = rib2Content.getString("@id");
1367 const int rib2ID = atoi(srib2ID.c_str());
1368 // const string rib2Name = rib2Content.getString("Name");
1369 const double length = rib2Content.getLength("Length");
1370 const double width = rib2Content.getLength("Width");
1371 const double thick = rib2Content.getLength("Thickness");
1372 const double width2 = rib2Content.getLength("Width2");
1373 const double thick2 = rib2Content.getLength("Thickness2");
1374 const double rotX = rib2Content.getLength("RotX");
1375 const double rotY = rib2Content.getLength("RotY");
1376 const double rotZ = rib2Content.getLength("RotZ");
1377 const double cX = rib2Content.getLength("PosX");
1378 const double cY = rib2Content.getLength("PosY");
1379 const double cZ = rib2Content.getLength("PosZ");
1380 const int number = atoi((rib2Content.getString("NDiv")).c_str());
1381
1382 const string solidName = "solidRib2" + to_string(rib2ID);
1383 const string logicalName = "logicalRib2" + to_string(rib2ID);
1384 G4Trd* trdShape = new G4Trd(solidName,
1385 0.5 * thick * CLHEP::cm,
1386 0.5 * thick2 * CLHEP::cm,
1387 0.5 * width * CLHEP::cm,
1388 0.5 * width2 * CLHEP::cm,
1389 0.5 * length * CLHEP::cm);
1390
1391 G4LogicalVolume* logicalV = new G4LogicalVolume(trdShape, medAluminum, logicalName, 0, 0, 0);
1392 if (rib2ID > 0)
1393 logicalV = new G4LogicalVolume(trdShape, medCopper, logicalName, 0, 0, 0);
1394
1395 logicalV->SetVisAttributes(m_VisAttributes.back());
1396
1397 const double phi = 360.0 / number;
1398
1399 G4RotationMatrix rot = G4RotationMatrix();
1400 G4ThreeVector arm(cX * CLHEP::cm, cY * CLHEP::cm, cZ * CLHEP::cm - thick * CLHEP::cm / 2.0);
1401
1402 rot.rotateX(rotX);
1403 rot.rotateY(rotY);
1404 rot.rotateZ(rotZ);
1405 for (int i = 0; i < number; ++i) {
1406 const string physicalName = "physicalRib2_" + to_string(rib2ID) + " " + to_string(i);
1407 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
1408 physicalName.c_str(), m_logicalCDC, false, rib2ID);
1409 rot.rotateZ(phi * CLHEP::deg);
1410 arm.rotateZ(phi * CLHEP::deg);
1411 }
1412
1413 }// rib2
1414
1415 const int nRib3s = content.getNumberNodes("Covers/Rib3");
1416 for (int iRib3 = 0; iRib3 < nRib3s; ++iRib3) {
1417 GearDir rib3Content(content);
1418 rib3Content.append((boost::format("/Covers/Rib3[%1%]/") % (iRib3 + 1)).str());
1419 const string srib3ID = rib3Content.getString("@id");
1420 const int rib3ID = atoi(srib3ID.c_str());
1421 // const string rib3Name = rib3Content.getString("Name");
1422 const double length = rib3Content.getLength("Length");
1423 const double width = rib3Content.getLength("Width");
1424 const double thick = rib3Content.getLength("Thickness");
1425 const double r = rib3Content.getLength("HoleR");
1426 const double cX = rib3Content.getLength("PosX");
1427 const double cY = rib3Content.getLength("PosY");
1428 const double cZ = rib3Content.getLength("PosZ");
1429 const double hX = rib3Content.getLength("HoleX");
1430 const double hY = rib3Content.getLength("HoleY");
1431 const double hZ = rib3Content.getLength("HoleZ");
1432 const int offset = atoi((rib3Content.getString("Offset")).c_str());
1433 const int number = atoi((rib3Content.getString("NDiv")).c_str());
1434
1435 const string logicalName = "logicalRib3" + to_string(rib3ID);
1436 G4VSolid* boxShape = new G4Box("Block",
1437 0.5 * length * CLHEP::cm,
1438 0.5 * width * CLHEP::cm,
1439 0.5 * thick * CLHEP::cm);
1440 G4VSolid* tubeShape = new G4Tubs("Hole",
1441 0.,
1442 r * CLHEP::cm,
1443 length * CLHEP::cm,
1444 0.,
1445 360. * CLHEP::deg);
1446 G4RotationMatrix rotsub = G4RotationMatrix();
1447 G4ThreeVector trnsub(cX * CLHEP::cm - hX * CLHEP::cm, cY * CLHEP::cm - hY * CLHEP::cm,
1448 cZ * CLHEP::cm - hZ * CLHEP::cm + 0.5 * thick * CLHEP::cm);
1449 G4VSolid* coolingBlock = new G4SubtractionSolid("Block-Hole",
1450 boxShape,
1451 tubeShape,
1452 G4Transform3D(rotsub,
1453 trnsub));
1454
1455 G4LogicalVolume* logicalV = new G4LogicalVolume(coolingBlock, medCopper, logicalName, 0, 0, 0);
1456
1457 logicalV->SetVisAttributes(m_VisAttributes.back());
1458
1459 const double phi = 360.0 / number;
1460
1461 G4RotationMatrix rot = G4RotationMatrix();
1462 G4ThreeVector arm(cX * CLHEP::cm, cY * CLHEP::cm, cZ * CLHEP::cm - thick * CLHEP::cm / 2.0);
1463
1464 if (offset) {
1465 rot.rotateZ(0.5 * phi * CLHEP::deg);
1466 arm.rotateZ(0.5 * phi * CLHEP::deg);
1467 }
1468 for (int i = 0; i < number; ++i) {
1469 const string physicalName = "physicalRib3_" + to_string(rib3ID) + " " + to_string(i);
1470 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
1471 physicalName.c_str(), m_logicalCDC, false, rib3ID);
1472 rot.rotateZ(phi * CLHEP::deg);
1473 arm.rotateZ(phi * CLHEP::deg);
1474 }
1475
1476 }// rib3
1477
1478 const int nRib4s = content.getNumberNodes("Covers/Rib4");
1479 for (int iRib4 = 0; iRib4 < nRib4s; ++iRib4) {
1480 GearDir rib4Content(content);
1481 rib4Content.append((boost::format("/Covers/Rib4[%1%]/") % (iRib4 + 1)).str());
1482 const string srib4ID = rib4Content.getString("@id");
1483 const int rib4ID = atoi(srib4ID.c_str());
1484 // const string rib4Name = rib4Content.getString("Name");
1485 const double length = rib4Content.getLength("Length");
1486 const double width = rib4Content.getLength("Width");
1487 const double thick = rib4Content.getLength("Thickness");
1488 const double length2 = rib4Content.getLength("Length2");
1489 const double width2 = rib4Content.getLength("Width2");
1490 const double thick2 = rib4Content.getLength("Thickness2");
1491 const double cX = rib4Content.getLength("PosX");
1492 const double cY = rib4Content.getLength("PosY");
1493 const double cZ = rib4Content.getLength("PosZ");
1494 const double hX = rib4Content.getLength("HoleX");
1495 const double hY = rib4Content.getLength("HoleY");
1496 const double hZ = rib4Content.getLength("HoleZ");
1497 const int offset = atoi((rib4Content.getString("Offset")).c_str());
1498 const int number = atoi((rib4Content.getString("NDiv")).c_str());
1499
1500 const string logicalName = "logicalRib4" + to_string(rib4ID);
1501 G4VSolid* baseShape = new G4Box("Base",
1502 0.5 * length * CLHEP::cm,
1503 0.5 * width * CLHEP::cm,
1504 0.5 * thick * CLHEP::cm);
1505 G4VSolid* sqShape = new G4Box("Sq",
1506 0.5 * length2 * CLHEP::cm,
1507 0.5 * width2 * CLHEP::cm,
1508 0.5 * thick2 * CLHEP::cm);
1509 G4RotationMatrix rotsub = G4RotationMatrix();
1510 double dzc = (hZ - thick2 / 2.) - (cZ - thick / 2.);
1511 G4ThreeVector trnsub(hX * CLHEP::cm - cX * CLHEP::cm,
1512 hY * CLHEP::cm - cY * CLHEP::cm,
1513 dzc * CLHEP::cm);
1514 G4VSolid* sqHoleBase = new G4SubtractionSolid("Base-Sq",
1515 baseShape,
1516 sqShape,
1517 G4Transform3D(rotsub,
1518 trnsub)
1519 );
1520
1521 G4LogicalVolume* logicalV = new G4LogicalVolume(sqHoleBase, medCopper, logicalName, 0, 0, 0);
1522 if (rib4ID < 19)
1523 logicalV = new G4LogicalVolume(sqHoleBase, medNEMA_G10_Plate, logicalName, 0, 0, 0);
1524
1525 logicalV->SetVisAttributes(m_VisAttributes.back());
1526
1527 const double phi = 360.0 / number;
1528
1529 G4RotationMatrix rot = G4RotationMatrix();
1530 G4ThreeVector arm(cX * CLHEP::cm, cY * CLHEP::cm, cZ * CLHEP::cm - thick * CLHEP::cm / 2.0);
1531
1532 if (offset) {
1533 rot.rotateZ(0.5 * phi * CLHEP::deg);
1534 arm.rotateZ(0.5 * phi * CLHEP::deg);
1535 }
1536 for (int i = 0; i < number; ++i) {
1537 const string physicalName = "physicalRib4_" + to_string(rib4ID) + " " + to_string(i);
1538 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
1539 physicalName.c_str(), m_logicalCDC, false, rib4ID);
1540 rot.rotateZ(phi * CLHEP::deg);
1541 arm.rotateZ(phi * CLHEP::deg);
1542 }
1543
1544 }// rib4
1545
1546 const int nRib5s = content.getNumberNodes("Covers/Rib5");
1547 for (int iRib5 = 0; iRib5 < nRib5s; ++iRib5) {
1548 GearDir rib5Content(content);
1549 rib5Content.append((boost::format("/Covers/Rib5[%1%]/") % (iRib5 + 1)).str());
1550 const string srib5ID = rib5Content.getString("@id");
1551 const int rib5ID = atoi(srib5ID.c_str());
1552 // const string rib5Name = rib5Content.getString("Name");
1553 const double dr = rib5Content.getLength("DeltaR");
1554 const double dz = rib5Content.getLength("DeltaZ");
1555 const double width = rib5Content.getLength("Width");
1556 const double thick = rib5Content.getLength("Thickness");
1557 const double rin = rib5Content.getLength("Rin");
1558 const double rotX = rib5Content.getLength("RotX");
1559 const double rotY = rib5Content.getLength("RotY");
1560 const double rotZ = rib5Content.getLength("RotZ");
1561 const double cX = rib5Content.getLength("PosX");
1562 const double cY = rib5Content.getLength("PosY");
1563 const double cZ = rib5Content.getLength("PosZ");
1564 const int offset = atoi((rib5Content.getString("Offset")).c_str());
1565 const int number = atoi((rib5Content.getString("NDiv")).c_str());
1566
1567 const string solidName = "solidRib5" + to_string(rib5ID);
1568 const string logicalName = "logicalRib5" + to_string(rib5ID);
1569 const double rmax = rin + thick;
1570 const double rmin = rin;
1571 const double dphi = 2. * atan2(dz, dr);
1572 const double ddphi = thick * tan(dphi) / rin;
1573 const double ddphi2 = width / 2. * width / 2. / (cX + dr) / rin;
1574 const double cphi = dphi - ddphi - ddphi2;
1575 G4Tubs* tubeShape = new G4Tubs(solidName,
1576 rmin * CLHEP::cm,
1577 rmax * CLHEP::cm,
1578 0.5 * width * CLHEP::cm,
1579 0.,
1580 cphi);
1581
1582 G4LogicalVolume* logicalV = new G4LogicalVolume(tubeShape, medAluminum, logicalName, 0, 0, 0);
1583
1584 logicalV->SetVisAttributes(m_VisAttributes.back());
1585
1586 const double phi = 360.0 / number;
1587
1588 G4RotationMatrix rot = G4RotationMatrix();
1589
1590 //G4ThreeVector arm(cX * CLHEP::cm, cY * CLHEP::cm, cZ * CLHEP::cm - thick * CLHEP::cm / 2.0);
1591 G4ThreeVector arm(cX * CLHEP::cm, cY * CLHEP::cm, cZ * CLHEP::cm - rin * CLHEP::cm - thick * CLHEP::cm);
1592
1593 rot.rotateX(rotX);
1594 rot.rotateY(rotY);
1595 rot.rotateZ(rotZ);
1596 if (offset) {
1597 rot.rotateZ(0.5 * phi * CLHEP::deg);
1598 arm.rotateZ(0.5 * phi * CLHEP::deg);
1599 }
1600 for (int i = 0; i < number; ++i) {
1601 const string physicalName = "physicalRib5_" + to_string(rib5ID) + " " + to_string(i);
1602 new G4PVPlacement(G4Transform3D(rot, arm), logicalV,
1603 physicalName.c_str(), m_logicalCDC, false, rib5ID);
1604 rot.rotateZ(phi * CLHEP::deg);
1605 arm.rotateZ(phi * CLHEP::deg);
1606 }
1607 }//rib5
1608
1609 }
1610
1611
1613 {
1614 G4Material* medAl = geometry::Materials::get("Al");
1615 G4double density = 1.000 * CLHEP::g / CLHEP::cm3;
1616 G4double a = 1.01 * CLHEP::g / CLHEP::mole;
1617 G4Element* elH = new G4Element("Hydrogen", "H", 1., a);
1618 a = 16.00 * CLHEP::g / CLHEP::mole;
1619 G4Element* elO = new G4Element("Oxygen", "O", 8., a);
1620 G4Material* medH2O = new G4Material("water", density, 2);
1621 medH2O->AddElement(elH, 2);
1622 medH2O->AddElement(elO, 1);
1623 G4Material* medCu = geometry::Materials::get("Cu");
1624 G4Material* medLV = geometry::Materials::get("CDCLVCable");
1625 G4Material* medFiber = geometry::Materials::get("CDCOpticalFiber");
1626 G4Material* medCAT7 = geometry::Materials::get("CDCCAT7");
1627 G4Material* medTRG = geometry::Materials::get("CDCOpticalFiberTRG");
1628 G4Material* medHV = geometry::Materials::get("CDCHVCable");
1629
1630 m_VisAttributes.push_back(new G4VisAttributes(true, G4Colour(0., 1., 0.)));
1631 for (const auto& cover : geom.getCovers()) {
1632 const int coverID = cover.getId();
1633 const string coverName = "cover" + to_string(coverID);
1634 const double rmin1 = cover.getRmin1();
1635 const double rmin2 = cover.getRmin2();
1636 const double rmax1 = cover.getRmax1();
1637 const double rmax2 = cover.getRmax2();
1638 const double thick = cover.getThick();
1639 const double posZ = cover.getZ();
1640
1641 /*if (coverID == 7 || coverID == 10) {
1642 createCone(rmin1, rmax1, rmin2, rmax2, thick, posZ, coverID, medAl, coverName);
1643 } else {
1644 createTube(rmin1, rmax1, thick, posZ, coverID, medAl, coverName);
1645 }*/
1646 // ID dependent material definition
1647 if (coverID < 23) {
1648 if (coverID == 7 || coverID == 10) {
1649 createCone(rmin1, rmax1, rmin2, rmax2, thick, posZ, coverID, medAl, coverName);
1650 } else {
1651 createTube(rmin1, rmax1, thick, posZ, coverID, medAl, coverName);
1652 }
1653 }
1654 if (coverID > 22 && coverID < 29)
1655 createTube(rmin1, rmax1, thick, posZ, coverID, medCu, coverName);
1656 if (coverID > 28 && coverID < 35)
1657 createTorus(rmin1, rmax1, thick, posZ, coverID, medCu, coverName);
1658 if (coverID > 34 && coverID < 41)
1659 createTorus(rmin1, rmax1, thick, posZ, coverID, medH2O, coverName);
1660 if (coverID == 45 || coverID == 46)
1661 createTube(rmin1, rmax1, thick, posZ, coverID, medLV, coverName);
1662 if (coverID == 47 || coverID == 48)
1663 createTube(rmin1, rmax1, thick, posZ, coverID, medFiber, coverName);
1664 if (coverID == 49 || coverID == 50)
1665 createTube(rmin1, rmax1, thick, posZ, coverID, medCAT7, coverName);
1666 if (coverID == 51 || coverID == 52)
1667 createTube(rmin1, rmax1, thick, posZ, coverID, medTRG, coverName);
1668 if (coverID == 53)
1669 createTube(rmin1, rmax1, thick, posZ, coverID, medHV, coverName);
1670 }
1671 }
1672
1674 {
1675 G4Material* medHV = geometry::Materials::get("CDCHVCable");
1676 G4Material* medFiber = geometry::Materials::get("CDCOpticalFiber");
1677 G4Material* medCAT7 = geometry::Materials::get("CDCCAT7");
1678 G4Material* medTRG = geometry::Materials::get("CDCOpticalFiberTRG");
1679
1680 m_VisAttributes.push_back(new G4VisAttributes(true, G4Colour(0., 1., 0.)));
1681 for (const auto& cover2 : geom.getCover2s()) {
1682 const int cover2ID = cover2.getId();
1683 const string cover2Name = "cover2" + to_string(cover2ID);
1684 const double rmin = cover2.getRmin();
1685 const double rmax = cover2.getRmax();
1686 const double phis = cover2.getPhis();
1687 const double dphi = cover2.getDphi();
1688 const double thick = cover2.getThick();
1689 const double posZ = cover2.getZ();
1690
1691 if (cover2ID < 11)
1692 createTube2(rmin, rmax, phis, dphi, thick, posZ, cover2ID, medHV, cover2Name);
1693 if (cover2ID > 10 && cover2ID < 14)
1694 createTube2(rmin, rmax, phis, dphi, thick, posZ, cover2ID, medFiber, cover2Name);
1695 if (cover2ID > 13 && cover2ID < 23)
1696 createTube2(rmin, rmax, phis, dphi, thick, posZ, cover2ID, medCAT7, cover2Name);
1697 if (cover2ID > 22 && cover2ID < 29)
1698 createTube2(rmin, rmax, phis, dphi, thick, posZ, cover2ID, medTRG, cover2Name);
1699 }
1700 }
1701
1702 void GeoCDCCreator::createCone(const double rmin1, const double rmax1,
1703 const double rmin2, const double rmax2,
1704 const double thick, const double posZ,
1705 const int id, G4Material* med,
1706 const string& name)
1707 {
1708 const string solidName = "solid" + name;
1709 const string logicalName = "logical" + name;
1710 const string physicalName = "physical" + name;
1711 G4Cons* coverConeShape = new G4Cons(solidName.c_str(), rmin1 * CLHEP::cm, rmax1 * CLHEP::cm,
1712 rmin2 * CLHEP::cm, rmax2 * CLHEP::cm, thick * CLHEP::cm / 2.0, 0.*CLHEP::deg, 360.*CLHEP::deg);
1713 G4LogicalVolume* coverCone = new G4LogicalVolume(coverConeShape, med,
1714 logicalName.c_str(), 0, 0, 0);
1715 coverCone->SetVisAttributes(m_VisAttributes.back());
1716 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, posZ * CLHEP::cm - thick * CLHEP::cm / 2.0), coverCone,
1717 physicalName.c_str(), m_logicalCDC, false, id);
1718
1719 }
1720
1721 void GeoCDCCreator::createTube(const double rmin, const double rmax,
1722 const double thick, const double posZ,
1723 const int id, G4Material* med,
1724 const string& name)
1725 {
1726 const string solidName = "solid" + name;
1727 const string logicalName = "logical" + name;
1728 const string physicalName = "physical" + name;
1729 G4Tubs* solidV = new G4Tubs(solidName.c_str(),
1730 rmin * CLHEP::cm,
1731 rmax * CLHEP::cm,
1732 thick * CLHEP::cm / 2.0,
1733 0.*CLHEP::deg,
1734 360.*CLHEP::deg);
1735 G4LogicalVolume* logicalV = new G4LogicalVolume(solidV, med,
1736 logicalName.c_str(), 0, 0, 0);
1737 logicalV->SetVisAttributes(m_VisAttributes.back());
1738 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, posZ * CLHEP::cm - thick * CLHEP::cm / 2.0), logicalV,
1739 physicalName.c_str(), m_logicalCDC, false, id);
1740
1741 }
1742
1743 void GeoCDCCreator::createBox(const double length, const double height,
1744 const double thick, const double x,
1745 const double y, const double z,
1746 const int id, G4Material* med,
1747 const string& name)
1748 {
1749 const string solidName = (boost::format("solid%1%%2%") % name % id).str();
1750 const string logicalName = (boost::format("logical%1%%2%") % name % id).str();
1751 const string physicalName = (boost::format("physical%1%%2%") % name % id).str();
1752 G4Box* boxShape = new G4Box(solidName.c_str(), 0.5 * length * CLHEP::cm,
1753 0.5 * height * CLHEP::cm,
1754 0.5 * thick * CLHEP::cm);
1755 G4LogicalVolume* logicalV = new G4LogicalVolume(boxShape, med,
1756 logicalName.c_str(), 0, 0, 0);
1757 logicalV->SetVisAttributes(m_VisAttributes.back());
1758 new G4PVPlacement(0, G4ThreeVector(x * CLHEP::cm, y * CLHEP::cm, z * CLHEP::cm - thick * CLHEP::cm / 2.0), logicalV,
1759 physicalName.c_str(), m_logicalCDC, false, id);
1760
1761 }
1762
1763 void GeoCDCCreator::createTorus(const double rmin1, const double rmax1,
1764 const double thick, const double posZ,
1765 const int id, G4Material* med,
1766 const string& name)
1767 {
1768 const string solidName = "solid" + name;
1769 const string logicalName = "logical" + name;
1770 const string physicalName = "physical" + name;
1771 const double rtor = (rmax1 + rmin1) / 2.;
1772 const double rmax = rmax1 - rtor;
1773 const double rmin = max((rmax - thick), 0.);
1774
1775 G4Torus* solidV = new G4Torus(solidName.c_str(),
1776 rmin * CLHEP::cm,
1777 rmax * CLHEP::cm,
1778 rtor * CLHEP::cm,
1779 0.*CLHEP::deg,
1780 360.*CLHEP::deg);
1781 G4LogicalVolume* logicalV = new G4LogicalVolume(solidV, med,
1782 logicalName.c_str(), 0, 0, 0);
1783 logicalV->SetVisAttributes(m_VisAttributes.back());
1784 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, posZ * CLHEP::cm), logicalV,
1785 physicalName.c_str(), m_logicalCDC, false, id);
1786
1787 }
1788
1789 void GeoCDCCreator::createTube2(const double rmin, const double rmax,
1790 const double phis, const double phie,
1791 const double thick, const double posZ,
1792 const int id, G4Material* med,
1793 const string& name)
1794 {
1795 const string solidName = "solid" + name;
1796 const string logicalName = "logical" + name;
1797 const string physicalName = "physical" + name;
1798 G4Tubs* solidV = new G4Tubs(solidName.c_str(),
1799 rmin * CLHEP::cm,
1800 rmax * CLHEP::cm,
1801 thick * CLHEP::cm / 2.0,
1802 phis,
1803 phie);
1804 G4LogicalVolume* logicalV = new G4LogicalVolume(solidV, med,
1805 logicalName.c_str(), 0, 0, 0);
1806 logicalV->SetVisAttributes(m_VisAttributes.back());
1807 new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, posZ * CLHEP::cm - thick * CLHEP::cm / 2.0), logicalV,
1808 physicalName.c_str(), m_logicalCDC, false, id);
1809
1810 }
1811
1812 void GeoCDCCreator::createMapper(G4LogicalVolume& topVolume)
1813 {
1815 if (!gcp.getMapperGeometry()) return;
1816
1817 const double xc = 0.5 * (-0.0002769 + 0.0370499) * CLHEP::cm;
1818 const double yc = 0.5 * (-0.0615404 + -0.108948) * CLHEP::cm;
1819 const double zc = 0.5 * (-35.3 + 48.5) * CLHEP::cm;
1820 //3 plates
1821 // const double plateWidth = 13.756 * CLHEP::cm;
1822 // const double plateThick = 1.203 * CLHEP::cm;
1823 // const double plateLength = 83.706 * CLHEP::cm;
1824 const double plateWidth = 13.8 * CLHEP::cm;
1825 const double plateThick = 1.2 * CLHEP::cm;
1826 const double plateLength = 83.8 * CLHEP::cm;
1827 const double phi = gcp.getMapperPhiAngle() * CLHEP::deg; //phi-angle in lab.
1828 // std::cout << "phi= " << phi << std::endl;
1829 // const double endRingRmin = 4.1135 * CLHEP::cm;
1830 // const double endRingRmax = 15.353 * CLHEP::cm;
1831 // const double endRingThick = 2.057 * CLHEP::cm;
1832 const double endPlateRmin = 4.0 * CLHEP::cm;
1833 const double endPlateRmax = 15.5 * CLHEP::cm;
1834 const double bwdEndPlateThick = 1.7 * CLHEP::cm;
1835 const double fwdEndPlateThick = 2.0 * CLHEP::cm;
1836
1837 G4Material* medAluminum = geometry::Materials::get("Al");
1838
1839 string name = "Plate";
1840 int pID = 0;
1841 G4Box* plateShape = new G4Box("solid" + name, .5 * plateWidth, .5 * plateThick, .5 * plateLength);
1842 G4LogicalVolume* logical0 = new G4LogicalVolume(plateShape, medAluminum, "logical" + name, 0, 0, 0);
1843 logical0->SetVisAttributes(m_VisAttributes.back());
1844 // const double x = .5 * plateWidth;
1845 const double x = xc + 0.5 * plateWidth;
1846 // const double y = endRingRmin;
1847 const double y = yc + endPlateRmin + 0.1 * CLHEP::cm;
1848 // double z = 2.871 * CLHEP::cm;
1849 G4ThreeVector xyz(x, y, zc);
1850 G4RotationMatrix rotM3 = G4RotationMatrix();
1851 xyz.rotateZ(phi);
1852 rotM3.rotateZ(phi);
1853 new G4PVPlacement(G4Transform3D(rotM3, xyz), logical0, "physical" + name, &topVolume, false, pID);
1854
1855 const double alf = 120. * CLHEP::deg;
1856 xyz.rotateZ(alf);
1857 rotM3.rotateZ(alf);
1858 new G4PVPlacement(G4Transform3D(rotM3, xyz), logical0, "physical" + name, &topVolume, false, pID + 1);
1859
1860 xyz.rotateZ(alf);
1861 rotM3.rotateZ(alf);
1862 new G4PVPlacement(G4Transform3D(rotM3, xyz), logical0, "physical" + name, &topVolume, false, pID + 2);
1863
1864 //Define 2 end-plates
1865 //bwd
1866 name = "BwdEndPlate";
1867 G4Tubs* BwdEndPlateShape = new G4Tubs("solid" + name, endPlateRmin, endPlateRmax, 0.5 * bwdEndPlateThick, 0., 360.*CLHEP::deg);
1868 G4LogicalVolume* logical1 = new G4LogicalVolume(BwdEndPlateShape, medAluminum, "logical" + name, 0, 0, 0);
1869 logical1->SetVisAttributes(m_VisAttributes.back());
1870 // z = -40.0105 * CLHEP::cm;
1871 double z = -35.3 * CLHEP::cm - 0.5 * bwdEndPlateThick;
1872 pID = 0;
1873 new G4PVPlacement(0, G4ThreeVector(xc, yc, z), logical1, "physical" + name, &topVolume, false, pID);
1874
1875 //fwd
1876 // z = 45.7525 * CLHEP::cm;
1877 // new G4PVPlacement(0, G4ThreeVector(0., 0., z), logical1, "physical" + name, &topVolume, false, pID + 1);
1878 name = "FwdEndPlate";
1879 G4Tubs* FwdEndPlateShape = new G4Tubs("solid" + name, endPlateRmin, endPlateRmax, 0.5 * fwdEndPlateThick, 0., 360.*CLHEP::deg);
1880 G4LogicalVolume* logical2 = new G4LogicalVolume(FwdEndPlateShape, medAluminum, "logical" + name, 0, 0, 0);
1881 logical2->SetVisAttributes(m_VisAttributes.back());
1882 z = 48.5 * CLHEP::cm + 0.5 * fwdEndPlateThick;
1883 new G4PVPlacement(0, G4ThreeVector(xc, yc, z), logical2, "physical" + name, &topVolume, false, pID);
1884 }
1885 }
1887}
DataType Z() const
access variable Z (= .at(2) without boundary check)
Definition B2Vector3.h:435
void SetX(DataType x)
set X/1st-coordinate
Definition B2Vector3.h:457
B2Vector3< DataType > Cross(const B2Vector3< DataType > &p) const
Cross product.
Definition B2Vector3.h:296
DataType X() const
access variable X (= .at(0) without boundary check)
Definition B2Vector3.h:431
DataType Y() const
access variable Y (= .at(1) without boundary check)
Definition B2Vector3.h:433
DataType Mag() const
The magnitude (rho in spherical coordinate system).
Definition B2Vector3.h:159
void SetY(DataType y)
set Y/2nd-coordinate
Definition B2Vector3.h:459
DataType Perp() const
The transverse component (R in cylindrical coordinate system).
Definition B2Vector3.h:200
The Class for BeamBackground Sensitive Detector.
int getNEndPlateLayers() const
Get the number of endplate layers.
double getR() const
Get Radius.
double getZbwd() const
Get bwd z-position.
std::vector< double > getRmin() const
Get the list of the Rmin coordinates.
The Class for CDC geometry.
Definition CDCGeometry.h:27
double getGlobalOffsetY() const
Get the global y offset of CDC wrt Belle2 coord.
FieldLayer getFieldLayer(int i) const
Get the i-th field layer.
OuterWall getOuterWall(int i) const
Get the i-th outer wall.
std::vector< Rib2 > getRib2s() const
Get the list of rib2s.
double getFiducialRmin() const
Get the fiducial Rmin of CDC sensitive volume.
std::vector< Rib > getRibs() const
Get the list of ribs.
double getGlobalOffsetX() const
Get the global x offset of CDC wrt Belle2 coord.
std::vector< InnerWall > getInnerWalls() const
Get the list of inner walls.
std::vector< Rib5 > getRib5s() const
Get the list of rib5s.
double getFiducialRmax() const
Get the fiducial Rmax of CDC sensitive volume.
int getNFieldWires() const
Get the number of field wires.
std::vector< Rib3 > getRib3s() const
Get the list of rib3s.
MotherVolume getMotherVolume() const
Get the mother volume geometry of CDC.
int getNSenseWires() const
Get the number of sense wires.
double getSenseDiameter() const
Get the diameter of sense wire.
std::vector< Frontend > getFrontends() const
Get the list of frontend layers.
InnerWall getInnerWall(int i) const
Get the i-th inner wall.
double getFieldDiameter() const
Get the diameter of field wire.
EndPlate getEndPlate(int i) const
Get the i-th endplate.
std::vector< Rib4 > getRib4s() const
Get the list of rib4s.
int getNSenseLayers() const
Get the number of sense layers.
std::vector< EndPlate > getEndPlates() const
Get the list of endplates.
double getFeedthroughLength() const
Get the length of feedthrough.
double getGlobalOffsetZ() const
Get the global z offset of CDC wrt Belle2 coord.
std::vector< OuterWall > getOuterWalls() const
Get the list of outer walls.
SenseLayer getSenseLayer(int i) const
Get i-th sense layer.
The Class for CDC Geometry Control Parameters.
bool getPrintMaterialTable() const
Get printMaterialTable flag.
double getMaterialDefinitionMode() const
Get material definition mode.
bool getMapperGeometry()
Get mapper geometry flag.
double getMapperPhiAngle()
Get mapper phi-angle.
static CDCGeoControlPar & getInstance()
Static method to get a reference to the CDCGeoControlPar instance.
The Class for CDC Geometry Parameters.
double fieldWireDiameter() const
Returns diameter of the field wire.
const B2Vector3D wireForwardPosition(uint layerId, int cellId, EWirePosition set=c_Base) const
Returns the forward position of the input sense wire.
double senseWireDiameter() const
Returns diameter of the sense wire.
const B2Vector3D wireBackwardPosition(uint layerId, int cellId, EWirePosition set=c_Base) const
Returns the backward position of the input sense wire.
unsigned nWiresInLayer(int layerId) const
Returns wire numbers in a layer.
static CDCGeometryPar & Instance(const CDCGeometry *=nullptr)
Static method to get a reference to the CDCGeometryPar instance.
The Class for CDC Sensitive Detector.
static CDCSimControlPar & getInstance()
Static method to get a reference to the CDCSimControlPar instance.
void createBox(const double length, const double height, const double thick, const double x, const double y, const double z, const int id, G4Material *med, const std::string &name)
Create G4Box.
~GeoCDCCreator()
The destructor of the GeoCDCCreator class.
std::vector< G4VisAttributes * > m_VisAttributes
Vector of pointers to G4VisAttributes.
std::vector< G4UserLimits * > m_userLimits
Vector of pointers to G4UserLimits.
BkgSensitiveDetector * m_bkgsensitive
Sensitive detector for background studies.
void createGeometry(const CDCGeometry &parameters, G4LogicalVolume &topVolume, geometry::GeometryTypes type)
Create G4 geometry of CDC.
void createTube2(const double rmin, const double rmax, const double phis, const double phie, const double thick, const double posZ, const int id, G4Material *med, const std::string &name)
Create G4Tube2.
void createTube(const double rmin, const double rmax, const double thick, const double posZ, const int id, G4Material *med, const std::string &name)
Create G4Tube.
GeoCDCCreator()
Constructor of the GeoCDCCreator class.
void createNeutronShields(const GearDir &content)
Create neutron shield from gearbox.
G4VPhysicalVolume * m_physicalCDC
CDC G4 physical volume.
void createTorus(const double rmin1, const double rmax1, const double thick, const double posZ, const int id, G4Material *med, const std::string &name)
Create G4Torus.
std::vector< BkgSensitiveDetector * > m_BkgSensitiveRib4
Sensitive detectors for background studies (rib4).
void createCone(const double rmin1, const double rmax1, const double rmin2, const double rmax2, const double thick, const double posz, const int id, G4Material *med, const std::string &name)
Create G4Cone.
void createCover2s(const GearDir &content)
Create CDC cover2s from gear box.
void createMapper(G4LogicalVolume &topVolume)
Create the B-field mapper geometry (tentative function)
G4LogicalVolume * m_logicalCDC
CDC G4 logical volume.
void createCovers(const GearDir &content)
Create CDC covers from gear box.
CDCSensitiveDetector * m_sensitive
Sensitive detector.
GearDir is the basic class used for accessing the parameter store.
Definition GearDir.h:31
void append(const std::string &path)
Append something to the current path, modifying the GearDir in place.
Definition GearDir.h:52
virtual std::string getString(const std::string &path="") const noexcept(false) override
Get the parameter path as a string.
Definition GearDir.h:69
double getLength(const std::string &path="") const noexcept(false)
Get the parameter path as a double converted to the standard length unit.
Definition Interface.h:259
static G4Material * get(const std::string &name)
Find given material.
Definition Materials.h:63
B2Vector3< double > B2Vector3D
typedef for common usage with double
Definition B2Vector3.h:516
double tan(double a)
tan for double
Definition beamHelpers.h:31
geometry::CreatorFactory< GeoCDCCreator > GeoCDCFactory("CDCCreator")
Register the GeoCreator.
Common code concerning the geometry representation of the detector.
Definition CreatorBase.h:25
GeometryTypes
Flag indicating the type of geometry to be used.
Abstract base class for different kinds of events.
STL namespace.
Very simple class to provide an easy way to register creators with the CreatorManager.