ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PHG4BbcDetector.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file PHG4BbcDetector.cc
1 #include "PHG4BbcDetector.h"
2 
3 #include "PHG4BbcDisplayAction.h"
4 
5 #include <phparameter/PHParameters.h>
6 
7 #include <g4main/PHG4Detector.h>
9 #include <g4main/PHG4Subsystem.h>
10 
11 #include <phool/recoConsts.h>
12 
13 #include <Geant4/G4LogicalVolume.hh>
14 #include <Geant4/G4Material.hh>
15 #include <Geant4/G4NistManager.hh>
16 #include <Geant4/G4PVPlacement.hh>
17 #include <Geant4/G4Polyhedra.hh>
18 #include <Geant4/G4RotationMatrix.hh> // for G4RotationMatrix
19 #include <Geant4/G4String.hh> // for G4String
20 #include <Geant4/G4SystemOfUnits.hh>
21 #include <Geant4/G4ThreeVector.hh>
22 #include <Geant4/G4Tubs.hh>
23 #include <Geant4/G4Types.hh> // for G4double, G4int
24 #include <Geant4/G4VPhysicalVolume.hh>
25 
26 #include <cmath>
27 #include <iostream> // for operator<<, endl, bas...
28 
29 class PHCompositeNode;
30 
31 PHG4BbcDetector::PHG4BbcDetector(PHG4Subsystem *subsys, PHCompositeNode *Node, PHParameters *params, const std::string &dnam)
32  : PHG4Detector(subsys, Node, dnam)
33  , m_DisplayAction(dynamic_cast<PHG4BbcDisplayAction *>(subsys->GetDisplayAction()))
34  , m_Params(params)
35  , m_ActiveFlag(m_Params->get_int_param("active"))
36  , m_SupportActiveFlag(m_Params->get_int_param("supportactive"))
37 {
38 }
39 
40 //_______________________________________________________________
41 //_______________________________________________________________
43 {
44  G4LogicalVolume *mylogvol = volume->GetLogicalVolume();
45 
46  if (m_ActiveFlag)
47  {
48  if (m_PhysLogicalVolSet.find(mylogvol) != m_PhysLogicalVolSet.end())
49  {
50  return 1;
51  }
52  }
54  {
55  if (m_SupportLogicalVolSet.find(mylogvol) != m_SupportLogicalVolSet.end())
56  {
57  return -2;
58  }
59  }
60  return 0;
61 }
62 
64 {
65  //std::cout << "PHG4BbcDetector::ConstructMe()" << std::endl;
66 
68 
69  //========================
70  // Build a single BBC PMT
71  //========================
72 
73  // BBC Detector Tube (BBCD) Logical Volume (contains all parts of PMT+Radiator)
75  G4double z_bbcd[] = {-6.15 * cm, 6.15 * cm};
76  G4double len_bbcd = z_bbcd[1] - z_bbcd[0];
77  G4double rin_bbcd[] = {0 * cm, 0 * cm};
78  G4double rout_bbcd[] = {1.4 * cm, 1.4 * cm};
79  G4Polyhedra *bbcd = new G4Polyhedra("bbcd", 0., 2 * M_PI, 6, 2, z_bbcd, rin_bbcd, rout_bbcd);
80  G4LogicalVolume *bbcd_lv = new G4LogicalVolume(bbcd, WorldMaterial, G4String("Bbc_tube"));
81 
82  //
83  // Place the BBCA, BBCQ, BBCP, BBCR and BBCH in to BBCD.
84  //
85 
86  // BBC Attachment Plate (BBCA)
87  const G4double z_bbca[] = {-0.5 * cm, -0.201 * cm, -0.2 * cm, 0.5 * cm};
88  const G4double rInner_bbca[] = {0.2 * cm, 0.2 * cm, 0.2 * cm, 0.2 * cm};
89  const G4double rOuter_bbca[] = {1.4 * cm, 1.4 * cm, 1.375 * cm, 1.375 * cm};
90  G4Polyhedra *bbca = new G4Polyhedra("bbca", 0., 2 * M_PI, 6, 4, z_bbca, rInner_bbca, rOuter_bbca);
91 
92  G4Material *Aluminum = GetDetectorMaterial("G4_Al");
93 
94  G4LogicalVolume *bbca_lv = new G4LogicalVolume(bbca, Aluminum, G4String("Bbc_attach_plate"));
95  m_SupportLogicalVolSet.insert(bbca_lv);
96  GetDisplayAction()->AddVolume(bbca_lv, "Bbc_attach_plate");
97 
98  G4double xpos = 0. * cm;
99  G4double ypos = 0. * cm;
100  G4double len_bbca = z_bbca[3] - z_bbca[0];
101  G4double zpos = z_bbcd[0] + len_bbca * 0.5;
102 
103  G4VPhysicalVolume *bbca_phys = new G4PVPlacement(nullptr, G4ThreeVector(xpos, ypos, zpos), bbca_lv, "BBCA", bbcd_lv, false, 0);
104  if (!bbca_phys) // more to prevent compiler warnings
105  {
106  std::cout << "placement of BBCA failed" << std::endl;
107  }
108 
109  G4Material *Quartz = GetDetectorMaterial("G4_SILICON_DIOXIDE");
110 
111  // BBC Quartz Radiator
112  const G4double z_bbcq[] = {-1.5 * cm, 1.5 * cm};
113  const G4double rInner_bbcq[] = {0., 0.};
114  const G4double rOuter_bbcq[] = {1.27 * cm, 1.27 * cm};
115  G4Polyhedra *bbcq = new G4Polyhedra("bbcq", 0., 2 * M_PI, 6, 2, z_bbcq, rInner_bbcq, rOuter_bbcq);
116 
117  G4LogicalVolume *bbcq_lv = new G4LogicalVolume(bbcq, Quartz, G4String("Bbc_quartz"));
118  GetDisplayAction()->AddVolume(bbcq_lv, "Bbc_quartz");
119 
120  xpos = 0. * cm;
121  ypos = 0. * cm;
122  G4double len_bbcq = z_bbcq[1] - z_bbcq[0];
123  zpos += len_bbca * 0.5 + len_bbcq * 0.5;
124 
125  G4VPhysicalVolume *bbcq_phys = new G4PVPlacement(nullptr, G4ThreeVector(xpos, ypos, zpos), bbcq_lv, "BBCQ", bbcd_lv, false, 0);
126  if (!bbcq_phys) // more to prevent compiler warnings
127  {
128  std::cout << "placement of BBCQ failed" << std::endl;
129  }
130 
131  // BBC PMT
132  const G4double len_bbcp = 4.4 * cm;
133  const G4double rInner_bbcp = 1.09 * cm;
134  const G4double rOuter_bbcp = 1.29 * cm;
135  G4Tubs *bbcp = new G4Tubs("bbcp", rInner_bbcp, rOuter_bbcp, len_bbcp * 0.5, 0 * deg, 360 * deg);
136 
137  G4LogicalVolume *bbcp_lv = new G4LogicalVolume(bbcp, Quartz, G4String("Bbc_PMT"));
138  GetDisplayAction()->AddVolume(bbcp_lv, "Bbc_PMT");
139 
140  xpos = 0. * cm;
141  ypos = 0. * cm;
142  zpos += len_bbcq * 0.5 + len_bbcp * 0.5;
143 
144  G4VPhysicalVolume *bbcp_phys = new G4PVPlacement(nullptr, G4ThreeVector(xpos, ypos, zpos), bbcp_lv, "BBCP", bbcd_lv, false, 0);
145  if (!bbcp_phys) // more to prevent compiler warnings
146  {
147  std::cout << "placement of BBCP failed" << std::endl;
148  }
149 
150  // BBC Breeder Module
151  const G4double len_bbcr = 3.9 * cm;
152  const G4double rInner_bbcr = 0.0 * cm;
153  const G4double rOuter_bbcr = 1.2 * cm;
154  G4Tubs *bbcr = new G4Tubs("bbcr", rInner_bbcr, rOuter_bbcr, len_bbcr * 0.5, 0 * deg, 360 * deg);
155 
156  G4int natoms;
157  G4int ncomponents;
158  G4double density;
159  G4Material *G10 = new G4Material("BBC_G10", density = 1.700 * g / cm3, ncomponents = 4);
161  G10->AddElement(G4NistManager::Instance()->FindOrBuildElement("Si"), natoms = 1);
162  G10->AddElement(G4NistManager::Instance()->FindOrBuildElement("O"), natoms = 2);
163  G10->AddElement(G4NistManager::Instance()->FindOrBuildElement("C"), natoms = 3);
164  G10->AddElement(G4NistManager::Instance()->FindOrBuildElement("H"), natoms = 3);
165 
166  G4LogicalVolume *bbcr_lv = new G4LogicalVolume(bbcr, G10, "Bbc_Breeder_Module");
167  GetDisplayAction()->AddVolume(bbcr_lv, "Bbc_Breeder_Module");
168 
169  xpos = 0. * cm;
170  ypos = 0. * cm;
171  zpos += len_bbcp * 0.5 + len_bbcr * 0.5;
172 
173  G4PVPlacement *plcmt = new G4PVPlacement(nullptr, G4ThreeVector(xpos, ypos, zpos), bbcr_lv, "BBCR", bbcd_lv, false, 0);
174  if (!plcmt) // more to prevent compiler warnings
175  {
176  std::cout << "placement of BBCR failed" << std::endl;
177  }
178 
179  // BBC Mu Metal Shield
180  const G4double z_bbch[] = {-5.65 * cm, 5.65 * cm};
181  const G4double rInner_bbch[] = {1.375 * cm, 1.375 * cm};
182  const G4double rOuter_bbch[] = {1.4 * cm, 1.4 * cm};
183  G4Polyhedra *bbch = new G4Polyhedra("bbch", 0., 2 * M_PI, 6, 2, z_bbch, rInner_bbch, rOuter_bbch);
184 
185  G4Material *MuMetal = manager->FindOrBuildMaterial("G4_STAINLESS-STEEL");
186 
187  G4LogicalVolume *bbch_lv = new G4LogicalVolume(bbch, MuMetal, G4String("Bbc_Shield"));
188 
189  xpos = 0. * cm;
190  ypos = 0. * cm;
191  G4double len_bbch = z_bbch[1] - z_bbch[0];
192  zpos = z_bbcd[0] + 0.3 * cm + len_bbch * 0.5;
193 
194  G4VPhysicalVolume *bbch_phys = new G4PVPlacement(nullptr, G4ThreeVector(xpos, ypos, zpos), bbch_lv, "BBCH", bbcd_lv, false, 0);
195  if (!bbch_phys) // more to prevent compiler warnings
196  {
197  std::cout << "placement of BBCH failed" << std::endl;
198  }
199 
200  // Locations of the 64 PMT tubes in an arm.
201  // These are the x,y for the south BBC.
202  // The north inverts the x coordinate (x -> -x)
203  // (NB: Should probably move this to a geometry object...)
204  float TubeLoc[64][2] = {
205  { -12.2976, 4.26 },
206  { -12.2976, 1.42 },
207  { -9.83805, 8.52 },
208  { -9.83805, 5.68 },
209  { -9.83805, 2.84 },
210  { -7.37854, 9.94 },
211  { -7.37854, 7.1 },
212  { -7.37854, 4.26 },
213  { -7.37854, 1.42 },
214  { -4.91902, 11.36 },
215  { -4.91902, 8.52 },
216  { -4.91902, 5.68 },
217  { -2.45951, 12.78 },
218  { -2.45951, 9.94 },
219  { -2.45951, 7.1 },
220  { 0, 11.36 },
221  { 0, 8.52 },
222  { 2.45951, 12.78 },
223  { 2.45951, 9.94 },
224  { 2.45951, 7.1 },
225  { 4.91902, 11.36 },
226  { 4.91902, 8.52 },
227  { 4.91902, 5.68 },
228  { 7.37854, 9.94 },
229  { 7.37854, 7.1 },
230  { 7.37854, 4.26 },
231  { 7.37854, 1.42 },
232  { 9.83805, 8.52 },
233  { 9.83805, 5.68 },
234  { 9.83805, 2.84 },
235  { 12.2976, 4.26 },
236  { 12.2976, 1.42 },
237  { 12.2976, -4.26 },
238  { 12.2976, -1.42 },
239  { 9.83805, -8.52 },
240  { 9.83805, -5.68 },
241  { 9.83805, -2.84 },
242  { 7.37854, -9.94 },
243  { 7.37854, -7.1 },
244  { 7.37854, -4.26 },
245  { 7.37854, -1.42 },
246  { 4.91902, -11.36 },
247  { 4.91902, -8.52 },
248  { 4.91902, -5.68 },
249  { 2.45951, -12.78 },
250  { 2.45951, -9.94 },
251  { 2.45951, -7.1 },
252  { 0, -11.36 },
253  { 0, -8.52 },
254  { -2.45951, -12.78 },
255  { -2.45951, -9.94 },
256  { -2.45951, -7.1 },
257  { -4.91902, -11.36 },
258  { -4.91902, -8.52 },
259  { -4.91902, -5.68 },
260  { -7.37854, -9.94 },
261  { -7.37854, -7.1 },
262  { -7.37854, -4.26 },
263  { -7.37854, -1.42 },
264  { -9.83805, -8.52 },
265  { -9.83805, -5.68 },
266  { -9.83805, -2.84 },
267  { -12.2976, -4.26 },
268  { -12.2976, -1.42 }
269  };
270 
271  //m_bbcz = m_Params->get_double_param("z") * cm;
272  m_bbcz = 250.0 * cm; // The front face of the quartz is at 250 cm
273 
274  const float tube_zpos = m_bbcz + len_bbcd / 2.0 - len_bbca;
275 
276  const int NPMT = 64; // No. PMTs per arm
277  G4RotationMatrix *arm_rot[2];
278  for (int iarm = 0; iarm < 2; iarm++)
279  {
280  arm_rot[iarm] = new G4RotationMatrix;
281  float xside = 1.0;
282  float zside = 1.0;
283  if (iarm == 0) // South Arm = 0
284  {
285  xside = 1.0;
286  zside = -1.0;
287  arm_rot[iarm]->rotateY(180 * deg);
288  }
289  else // North Arm = 1
290  {
291  xside = -1.0;
292  zside = 1.0;
293  }
294 
295  // Add BBC PMT's
296  for (int itube = 0; itube < NPMT; itube++)
297  {
298  // Full PMT Housing with Active Quartz Cerenkov Radiators
299  float tube_xpos = xside * TubeLoc[itube][0] * cm;
300  float tube_ypos = TubeLoc[itube][1] * cm;
301  new G4PVPlacement(arm_rot[iarm], G4ThreeVector( tube_xpos, tube_ypos, zside*tube_zpos ),
302  bbcd_lv, "BBCD", logicWorld, false, iarm * NPMT + itube, OverlapCheck());
303  }
304  }
305 
306  // Now Build the BBC Housing
307 
308  // BBC Outer and Inner Cylindrical Shells
309  G4Tubs *bbc_outer_shell = new G4Tubs("bbc_outer_shell", 14.9 * cm, 15 * cm, 11.5 * cm, 0, 2 * M_PI);
310  G4LogicalVolume *bbc_outer_shell_lv = new G4LogicalVolume(bbc_outer_shell, Aluminum, G4String("Bbc_Outer_Shell"));
311  GetDisplayAction()->AddVolume(bbc_outer_shell_lv, "Bbc_Outer_Shell");
312 
313  G4Tubs *bbc_inner_shell = new G4Tubs("bbc_inner_shell", 5.0 * cm, 5.5 * cm, 11.5 * cm, 0, 2 * M_PI);
314  G4LogicalVolume *bbc_inner_shell_lv = new G4LogicalVolume(bbc_inner_shell, Aluminum, G4String("Bbc_Inner_Shell"));
315  GetDisplayAction()->AddVolume(bbc_inner_shell_lv, "Bbc_Inner_Shell");
316 
317  G4VPhysicalVolume *outer_shell_vol[2] = {0};
318  G4VPhysicalVolume *inner_shell_vol[2] = {0};
319 
320  // Place South Shells
321  outer_shell_vol[0] = new G4PVPlacement(nullptr, G4ThreeVector(0, 0, (-250 + 1.0 - 11.5) * cm),
322  bbc_outer_shell_lv, "BBC_OUTER_SHELL", logicWorld, false, 0, OverlapCheck());
323  inner_shell_vol[0] = new G4PVPlacement(nullptr, G4ThreeVector(0, 0, (-250 + 1.0 - 11.5) * cm),
324  bbc_inner_shell_lv, "BBC_INNER_SHELL", logicWorld, false, 0, OverlapCheck());
325 
326  // Place North Shells
327  outer_shell_vol[1] = new G4PVPlacement(nullptr, G4ThreeVector(0, 0, (250 - 1.0 + 11.5) * cm),
328  bbc_outer_shell_lv, "BBC_OUTER_SHELL", logicWorld, false, 1, OverlapCheck());
329  inner_shell_vol[1] = new G4PVPlacement(nullptr, G4ThreeVector(0, 0, (250 - 1.0 + 11.5) * cm),
330  bbc_inner_shell_lv, "BBC_INNER_SHELL", logicWorld, false, 0, OverlapCheck());
331  // this is more to prevent compiler warnings about unused variables
332  if (!outer_shell_vol[0] || !outer_shell_vol[1] || !inner_shell_vol[0] || !inner_shell_vol[1])
333  {
334  std::cout << "problem placing BBC Sheels" << std::endl;
335  }
336 
337  // BBC Front and Back Plates
338  G4Tubs *bbc_plate = new G4Tubs("bbc_fplate", 5 * cm, 15 * cm, 0.5 * cm, 0, 2 * M_PI);
339  G4LogicalVolume *bbc_plate_lv = new G4LogicalVolume(bbc_plate, Aluminum, G4String("Bbc_Cover_Plates"));
340  GetDisplayAction()->AddVolume(bbc_plate_lv, "Bbc_Cover_Plates");
341 
342  G4VPhysicalVolume *fplate_vol[2] = {0}; // Front Plates
343  G4VPhysicalVolume *bplate_vol[2] = {0}; // Back Plates
344 
345  // Place South Plates
346  fplate_vol[0] = new G4PVPlacement(nullptr, G4ThreeVector(0, 0, (-250 + 1.0 + 0.5) * cm),
347  bbc_plate_lv, "BBC_FPLATE", logicWorld, false, 0, OverlapCheck());
348  bplate_vol[0] = new G4PVPlacement(nullptr, G4ThreeVector(0, 0, (-250 + 1.0 + 0.5 - 24.0) * cm),
349  bbc_plate_lv, "BBC_BPLATE", logicWorld, false, 0, OverlapCheck());
350 
351  // Place North Plates
352  fplate_vol[1] = new G4PVPlacement(nullptr, G4ThreeVector(0, 0, (250 - 1.0 - 0.5) * cm),
353  bbc_plate_lv, "BBC_FPLATE", logicWorld, false, 1, OverlapCheck());
354  bplate_vol[1] = new G4PVPlacement(nullptr, G4ThreeVector(0, 0, (250 - 1.0 - 0.5 + 24.0) * cm),
355  bbc_plate_lv, "BBC_BPLATE", logicWorld, false, 0, OverlapCheck());
356 
357  // Place BBC Cables
358  G4Material* Cu = manager->FindOrBuildMaterial("G4_Cu");
359  const G4double len_cable = 120*cm;
360  const G4double r_CableConductor = 0.09525*cm;
361  G4Tubs *bbc_cablecond = new G4Tubs("bbc_cablecond",0.,r_CableConductor,len_cable*0.5,0*deg,360*deg);
362 
363  G4LogicalVolume *bbc_cablecond_lv = new G4LogicalVolume(bbc_cablecond, Cu, G4String("Bbc_CableCond"));
364 
365  const G4double rIn_CableShield = 0.302876*cm;
366  const G4double rOut_CableShield = 0.3175*cm;
367  G4Tubs *bbc_cableshield = new G4Tubs("bbc_cableshield",rIn_CableShield,rOut_CableShield,len_cable*0.5,0*deg,360*deg);
368 
369  G4LogicalVolume *bbc_cableshield_lv = new G4LogicalVolume(bbc_cableshield, Cu, G4String("Bbc_CableShield"));
370 
371  ypos = len_cable/2 + 5*cm;
372 
373  // For now we make this vertical, but they should be slanted toward the endcap
374  G4RotationMatrix *rot_cable = new G4RotationMatrix();
375  rot_cable->rotateX( 90*deg );
376 
377  int icable = 0;
378 
379  for (int iarm=0; iarm<2; iarm++)
380  {
381  float zsign = -1.0;
382  if ( iarm == 1 )
383  {
384  zsign = 1.0;
385  }
386 
387  for (int iring=1; iring<5; iring++)
388  {
389  float ring_radius = iring*0.67*cm;
390  int ncables = 2*M_PI*ring_radius/(0.635*cm);
391  double dphi = 2*M_PI/ncables;
392 
393  //G4cout << "BBC_CABLE " << iring << "\t" << ring_radius << "\t" << ncables << "\t" << dphi*180/3.14 << G4endl;
394 
395  // place cables in ring
396  for (int ic=0; ic<ncables; ic++)
397  {
398  xpos = cos(dphi*ic)*ring_radius;
399  zpos = sin(dphi*ic)*ring_radius + zsign*(m_bbcz + 30*cm);
400  // Place Inner Conductor
401  new G4PVPlacement(rot_cable, G4ThreeVector(xpos, ypos, zpos), bbc_cablecond_lv, "BBC_Cable_Cond", logicWorld, false, icable, OverlapCheck());
402  // Place Shield
403  new G4PVPlacement(rot_cable, G4ThreeVector(xpos, ypos, zpos), bbc_cableshield_lv, "BBC_Cable_Shield", logicWorld, false, icable++, OverlapCheck());
404  }
405  }
406  }
407 
408  // this is more to prevent compiler warnings about unused variables
409  if ( !fplate_vol[0] || !fplate_vol[1] || !bplate_vol[0] || !bplate_vol[1] )
410  {
411  std::cout << "problem placing BBC Sheets" << std::endl;
412  }
413 
414  // bbcq is the active detector element
415  m_PhysLogicalVolSet.insert(bbcq_lv);
416 
417  // GetDisplayAction()->AddVolume(bbc_plate_lv, "Bbc_plate");
418  // GetDisplayAction()->AddVolume(bbc_outer_shell_lv, "Bbc_oshell");
419  // GetDisplayAction()->AddVolume(bbc_inner_shell_lv, "Bbc_ishell");
420 
421  std::cout << "INSIDE BBC" << std::endl;
422 
423  return;
424 }
425 
426 void PHG4BbcDetector::Print(const std::string &what) const
427 {
428  std::cout << "Bbc Detector:" << std::endl;
429  if (what == "ALL" || what == "VOLUME")
430  {
431  std::cout << "Version 0.1" << std::endl;
432  }
433  return;
434 }