Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PHG4InttDetector.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file PHG4InttDetector.cc
1 #include "PHG4InttDetector.h"
2 
3 #include "PHG4InttDefs.h" // for SEGMENTATION_Z
6 
8 
10 
11 #include <phparameter/PHParameters.h>
12 #include <phparameter/PHParametersContainer.h>
13 
14 #include <g4main/PHG4Detector.h> // for PHG4Detector
15 #include <g4main/PHG4DisplayAction.h> // for PHG4DisplayAction
16 #include <g4main/PHG4Subsystem.h> // for PHG4Subsystem
17 
18 #include <phool/PHCompositeNode.h>
19 #include <phool/PHIODataNode.h>
20 #include <phool/PHNode.h> // for PHNode
21 #include <phool/PHNodeIterator.h> // for PHNodeIterator
22 #include <phool/PHObject.h> // for PHObject
23 #include <phool/getClass.h>
24 #include <phool/phool.h> // for PHWHERE
25 #include <phool/recoConsts.h>
26 
27 #include <TSystem.h>
28 
29 #include <Geant4/G4Box.hh>
30 #include <Geant4/G4GenericTrap.hh>
31 #include <Geant4/G4LogicalVolume.hh>
32 #include <Geant4/G4PVParameterised.hh>
33 #include <Geant4/G4PVPlacement.hh>
34 #include <Geant4/G4RotationMatrix.hh> // for G4RotationMatrix
35 #include <Geant4/G4String.hh> // for G4String
36 #include <Geant4/G4SubtractionSolid.hh>
37 #include <Geant4/G4SystemOfUnits.hh>
38 #include <Geant4/G4ThreeVector.hh> // for G4ThreeVector
39 #include <Geant4/G4Transform3D.hh> // for G4Transform3D
40 #include <Geant4/G4Tubs.hh>
41 #include <Geant4/G4TwoVector.hh> // for G4TwoVector
42 #include <Geant4/G4VPhysicalVolume.hh> // for G4VPhysicalVolume
43 #include <Geant4/geomdefs.hh> // for kZAxis
44 
45 #pragma GCC diagnostic push
46 #pragma GCC diagnostic ignored "-Wshadow"
47 #include <boost/format.hpp>
48 #pragma GCC diagnostic pop
49 
50 #include <algorithm> // for fill_n
51 #include <array>
52 #include <cassert> // for assert
53 #include <cmath>
54 #include <cstdlib> // for exit, NULL
55 #include <iostream> // for operator<<, basic...
56 
57 class G4VPVParameterisation;
58 class G4VSolid;
59 
60 PHG4InttDetector::PHG4InttDetector(PHG4Subsystem *subsys, PHCompositeNode *Node, PHParametersContainer *parameters, const std::string &dnam, const std::pair<std::vector<std::pair<int, int>>::const_iterator, std::vector<std::pair<int, int>>::const_iterator> &layer_b_e)
61  : PHG4Detector(subsys, Node, dnam)
62  , m_DisplayAction(dynamic_cast<PHG4InttDisplayAction *>(subsys->GetDisplayAction()))
63  , m_ParamsContainer(parameters)
64  , m_LayerBeginEndIteratorPair(layer_b_e)
65 {
66  for (auto layeriter = m_LayerBeginEndIteratorPair.first; layeriter != m_LayerBeginEndIteratorPair.second; ++layeriter)
67  {
68  int layer = layeriter->second;
69  const PHParameters *par = m_ParamsContainer->GetParameters(layer);
70  m_IsActiveMap.insert(std::make_pair(layer, par->get_int_param("active")));
71  m_IsAbsorberActiveMap.insert(std::make_pair(layer, par->get_int_param("absorberactive")));
72  }
74  m_IsSupportActive = par->get_int_param("supportactive");
75  m_IsEndcapActive = par->get_int_param("endcap_ring_enabled");
76  std::fill_n(&m_PosZ[0][0], sizeof(m_PosZ) / sizeof(double), NAN);
77  std::fill_n(m_SensorRadius, sizeof(m_SensorRadius) / sizeof(double), NAN);
78  std::fill_n(m_StripOffsetX, sizeof(m_StripOffsetX) / sizeof(double), NAN);
79 }
80 
81 //_______________________________________________________________
83 {
84  // Is this volume one of the sensor strips?
85  // just checking if the pointer to the logical volume is in the set
86  // of our active/active absorber ones makes sure we are in an active volume
87  // name parsing is a bad idea since this is called for all steps
88  // and we would have to trust that people give different names
89  // to their volumes
90  G4LogicalVolume *logvol = volume->GetLogicalVolume();
91  if (!m_PassiveVolumeTuple.empty() && m_PassiveVolumeTuple.find(logvol) != m_PassiveVolumeTuple.end())
92  {
93  return -1;
94  }
95  if (m_ActiveLogVols.find(logvol) != m_ActiveLogVols.end())
96  {
97  return 1;
98  }
99 
100  return 0;
101 }
102 
103 void PHG4InttDetector::ConstructMe(G4LogicalVolume *logicWorld)
104 {
105  if (Verbosity() > 0)
106  {
107  std::cout << "PHG4InttDetector::Construct called for layers " << std::endl;
108  for (auto layeriter = m_LayerBeginEndIteratorPair.first; layeriter != m_LayerBeginEndIteratorPair.second; ++layeriter)
109  {
110  std::cout << "layer " << layeriter->second << std::endl;
111  }
112  }
113  // the tracking layers are placed directly in the world volume, since some layers are (touching) double layers
114  ConstructIntt(logicWorld);
115 
116  // This object provides the strip center locations when given the ladder segment and strip internal cordinates in the sensor
117  AddGeometryNode();
118  return;
119 }
120 
121 int PHG4InttDetector::ConstructIntt(G4LogicalVolume *trackerenvelope)
122 {
123  recoConsts *rc = recoConsts::instance(); // use for worldmaterial in a few places
124  // We have an arbitray number of layers (nlayer_) up to 8
125  // We have 2 types of ladders (vertical strips and horizontal strips)
126  // We have 2 types of sensors (inner and outer)
127  std::array<std::array<double, 2>, 8> hdi_z_arr;
128  // we loop over layers. All layers have only one laddertype
129  for (auto layeriter = m_LayerBeginEndIteratorPair.first; layeriter != m_LayerBeginEndIteratorPair.second; ++layeriter)
130  {
131  int inttlayer = layeriter->second;
132  // get the parameters for this layer
133  const PHParameters *params1 = m_ParamsContainer->GetParameters(inttlayer);
134  const int laddertype = params1->get_int_param("laddertype");
135  const double offsetphi = (params1->get_double_param("offsetphi") * deg) / rad; // use rad internally
136  double offsetrot = (params1->get_double_param("offsetrot") * deg) / rad; // offsetrot is specified in deg, we convert to rad here
137  m_SensorRadius[inttlayer] = params1->get_double_param("sensor_radius") * cm;
138  const int nladders_layer = params1->get_int_param("nladder");
139 
140  // Look up all remaining parameters by the laddertype for this layer
141  const PHParameters *params = m_ParamsContainer->GetParameters(laddertype);
142  const double strip_x = params->get_double_param("strip_x") * cm;
143  const double strip_y = params->get_double_param("strip_y") * cm;
144  const int nstrips_phi_sensor = params->get_int_param("nstrips_phi_sensor");
145  const double sensor_offset_y = params->get_double_param("sensor_offset_y") * cm;
146  const double hdi_y = params->get_double_param("hdi_y") * cm;
147  double hdi_kapton_x = params->get_double_param("hdi_kapton_x") * cm;
148  double hdi_copper_x = params->get_double_param("hdi_copper_x") * cm;
149  double fphx_x = params->get_double_param("fphx_x") * cm;
150  double fphx_y = params->get_double_param("fphx_y") * cm;
151  double fphx_z = params->get_double_param("fphx_z") * cm;
152  double fphx_offset_z = params->get_double_param("fphx_offset_z") * cm;
153 
154  double si_glue_x = params->get_double_param("si_glue_x") * cm;
155  double fphx_glue_x = params->get_double_param("fphx_glue_x") * cm;
156  double halfladder_inside_z = params->get_double_param("halfladder_inside_z") * cm;
157 
158  if (Verbosity() > 0)
159  {
160  std::cout << "Constructing Intt layer: " << std::endl;
161  std::cout << " layer " << inttlayer << " laddertype " << laddertype << " nladders_layer " << nladders_layer
162  << " sensor_radius " << m_SensorRadius[inttlayer] << " offsetphi " << offsetphi << " rad "
163  << " offsetphi " << offsetphi * rad / deg << " deg "
164  << std::endl;
165  }
166  // We loop over inner, then outer (wrt the beam-axis), sensors, where itype specifies the inner or outer sensor
167  // The rest of this loop will construct and put in place a section of a ladder corresponding to the Z range of this sensor only
168  for (int itype = 0; itype < 2; ++itype)
169  {
170  double strip_z;
171  int nstrips_z_sensor;
172  switch (itype)
173  {
174  case 0:
175  strip_z = params->get_double_param("strip_z_0") * cm;
176  nstrips_z_sensor = params->get_int_param("nstrips_z_sensor_0");
177  break;
178  case 1:
179  strip_z = params->get_double_param("strip_z_1") * cm;
180  nstrips_z_sensor = params->get_int_param("nstrips_z_sensor_1");
181  break;
182  default:
183  std::cout << "invalid itype " << itype << std::endl;
184  exit(1);
185  }
186 
187  // ----- Step 1 ------------------------------------------------------
188  // We make the volumes for Si-sensor, FPHX, HDI, and stave components
189  // We add them to the ladder later
190  //====================================================================
191 
192  // Create Si-sensor active volume
193  const double siactive_x = strip_x;
194  const double siactive_y = strip_y * nstrips_phi_sensor;
195  const double siactive_z = strip_z * nstrips_z_sensor;
196  G4VSolid *siactive_box = new G4Box((boost::format("siactive_box_%d_%d") % inttlayer % itype).str(), siactive_x / 2, siactive_y / 2., siactive_z / 2.);
197  G4LogicalVolume *siactive_volume = new G4LogicalVolume(siactive_box, GetDetectorMaterial("G4_Si"),
198  boost::str(boost::format("siactive_volume_%d_%d") % inttlayer % itype).c_str(), 0, 0, 0);
199  if ((m_IsActiveMap.find(inttlayer))->second > 0)
200  {
201  m_ActiveLogVols.insert(siactive_volume);
202  }
203  m_DisplayAction->AddVolume(siactive_volume, "SiActive");
204  // We do not subdivide the sensor in G4. We will assign hits to strips in the stepping action, using the geometry object
205 
206  // Si-sensor full (active+inactive) area
207  const double sifull_x = siactive_x;
208  const double sifull_y = siactive_y + 2.0 * params->get_double_param("sensor_edge_phi") * cm;
209  const double sifull_z = siactive_z + 2.0 * params->get_double_param("sensor_edge_z") * cm;
210  G4VSolid *sifull_box = new G4Box((boost::format("sifull_box_%d_%d") % inttlayer % itype).str(), sifull_x / 2., sifull_y / 2.0, sifull_z / 2.0);
211 
212  // Si-sensor inactive area
213  G4VSolid *siinactive_box = new G4SubtractionSolid((boost::format("siinactive_box_%d_%d") % inttlayer % itype).str(),
214  sifull_box, siactive_box, 0, G4ThreeVector(0, 0, 0));
215  G4LogicalVolume *siinactive_volume = new G4LogicalVolume(siinactive_box, GetDetectorMaterial("G4_Si"),
216  (boost::format("siinactive_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
217 
218  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
219  {
220  m_PassiveVolumeTuple.insert(std::make_pair(siinactive_volume, std::make_tuple(inttlayer, PHG4InttDefs::SI_INACTIVE)));
221  }
222  m_DisplayAction->AddVolume(siinactive_volume, "SiInActive");
223 
224  // Glue for Si-sensor full area
225  G4VSolid *si_glue_box = new G4Box((boost::format("si_glue_box_%d_%d") % inttlayer % itype).str(), si_glue_x / 2., sifull_y / 2.0, sifull_z / 2.0);
226 
227  G4LogicalVolume *si_glue_volume = new G4LogicalVolume(si_glue_box, GetDetectorMaterial("SilverEpoxyGlue_INTT"),
228  (boost::format("si_glue_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
229 
230  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
231  {
232  m_PassiveVolumeTuple.insert(std::make_pair(siinactive_volume, std::make_tuple(inttlayer, PHG4InttDefs::SI_GLUE)));
233  }
234  m_DisplayAction->AddVolume(si_glue_volume, "SiGlue");
235 
236  // Make the HDI Kapton and copper volumes
237  // This makes HDI volumes that matche this sensor in Z length
238  const double hdi_z = sifull_z + params->get_double_param("hdi_edge_z") * cm;
239  hdi_z_arr[inttlayer][itype] = hdi_z;
240  G4VSolid *hdi_kapton_box = new G4Box((boost::format("hdi_kapton_box_%d_%d") % inttlayer % itype).str(), hdi_kapton_x / 2., hdi_y / 2., hdi_z / 2.0);
241  G4LogicalVolume *hdi_kapton_volume = new G4LogicalVolume(hdi_kapton_box, GetDetectorMaterial("G4_KAPTON"),
242  (boost::format("hdi_kapton_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
243 
244  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
245  {
246  m_PassiveVolumeTuple.insert(std::make_pair(hdi_kapton_volume, std::make_tuple(inttlayer, PHG4InttDefs::HDI_KAPTON)));
247  }
248  G4VSolid *hdi_copper_box = new G4Box((boost::format("hdi_copper_box_%d_%d") % inttlayer % itype).str(), hdi_copper_x / 2., hdi_y / 2., hdi_z / 2.0);
249  G4LogicalVolume *hdi_copper_volume = new G4LogicalVolume(hdi_copper_box, GetDetectorMaterial("G4_Cu"),
250  (boost::format("hdi_copper_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
251  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
252  {
253  m_PassiveVolumeTuple.insert(std::make_pair(hdi_copper_volume, std::make_tuple(inttlayer, PHG4InttDefs::HDI_COPPER)));
254  }
255  m_DisplayAction->AddVolume(hdi_kapton_volume, "HdiKapton");
256  m_DisplayAction->AddVolume(hdi_copper_volume, "HdiCopper");
257 
258  // This is the part of the HDI that extends beyond the sensor inside the endcap ring
259  const double hdiext_z = (itype == 0) ? 0.000001 : halfladder_inside_z - hdi_z_arr[inttlayer][0] - hdi_z; // need to assign nonzero value for itype=0
260  G4VSolid *hdiext_kapton_box = new G4Box((boost::format("hdiext_kapton_box_%d_%s") % inttlayer % itype).str(),
261  hdi_kapton_x / 2., hdi_y / 2., hdiext_z / 2.0);
262  G4LogicalVolume *hdiext_kapton_volume = new G4LogicalVolume(hdiext_kapton_box, GetDetectorMaterial("G4_KAPTON"), // was "FPC"
263  (boost::format("hdiext_kapton_%d_%s") % inttlayer % itype).str(), 0, 0, 0);
264  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
265  {
266  m_PassiveVolumeTuple.insert(std::make_pair(hdiext_kapton_volume, std::make_tuple(inttlayer, PHG4InttDefs::HDIEXT_KAPTON)));
267  }
268  G4VSolid *hdiext_copper_box = new G4Box((boost::format("hdiext_copper_box_%d_%s") % inttlayer % itype).str(),
269  hdi_copper_x / 2., hdi_y / 2., hdiext_z / 2.0);
270  G4LogicalVolume *hdiext_copper_volume = new G4LogicalVolume(hdiext_copper_box, GetDetectorMaterial("G4_Cu"),
271  (boost::format("hdiext_copper_%d_%s") % inttlayer % itype).str(), 0, 0, 0);
272  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
273  {
274  m_PassiveVolumeTuple.insert(std::make_pair(hdiext_copper_volume, std::make_tuple(inttlayer, PHG4InttDefs::HDIEXT_COPPER)));
275  }
276  m_DisplayAction->AddVolume(hdiext_kapton_volume, "HdiKapton");
277  m_DisplayAction->AddVolume(hdiext_copper_volume, "HdiCopper");
278 
279  // FPHX
280  G4VSolid *fphx_box = new G4Box((boost::format("fphx_box_%d_%d") % inttlayer % itype).str(), fphx_x / 2., fphx_y / 2., fphx_z / 2.);
281  G4LogicalVolume *fphx_volume = new G4LogicalVolume(fphx_box, GetDetectorMaterial("G4_Si"),
282  (boost::format("fphx_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
283  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
284  {
285  m_PassiveVolumeTuple.insert(std::make_pair(fphx_volume, std::make_tuple(inttlayer, PHG4InttDefs::FPHX)));
286  }
287  m_DisplayAction->AddVolume(fphx_volume, "FPHX");
288 
289  const double gap_sensor_fphx = params->get_double_param("gap_sensor_fphx") * cm;
290 
291  // FPHX Container
292  // make a container for the FPHX chips needed for this sensor, and then place them in the container
293  G4VSolid *fphxcontainer_box = new G4Box((boost::format("fphxcontainer_box_%d_%d") % inttlayer % itype).str(),
294  fphx_x / 2., fphx_y / 2., hdi_z / 2.);
295  G4LogicalVolume *fphxcontainer_volume = new G4LogicalVolume(fphxcontainer_box, GetDetectorMaterial(rc->get_StringFlag("WorldMaterial")),
296  (boost::format("fphxcontainer_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
297  m_DisplayAction->AddVolume(fphxcontainer_volume, "FPHXContainer");
298 
299  // Install multiple FPHX volumes in the FPHX container volume
300  // one FPHX chip per cell - each cell is 128 channels
301  const double fphx_offsetx = 0.;
302  const double fphx_offsety = 0.;
303  int ncopy;
304  double offsetz, cell_length_z;
305 
306  if (laddertype == PHG4InttDefs::SEGMENTATION_Z) // vertical strips
307  {
308  // For laddertype 0, we have 5 cells per sensor, but the strips are vertical, so we have to treat it specially
309  ncopy = nstrips_z_sensor / 128.0;
310  }
311  else if (laddertype == PHG4InttDefs::SEGMENTATION_PHI)
312  {
313  ncopy = nstrips_z_sensor;
314  }
315  else
316  {
317  std::cout << PHWHERE << "invalid laddertype " << laddertype << std::endl;
318  gSystem->Exit(1);
319  // this is just to make the optimizer happy which otherwise complains about possibly
320  // uninitialized variables. It doesn't know gSystem->Exit(1) quits,
321  // this exit here terminates the program for it
322  exit(1);
323  }
324  cell_length_z = strip_z * nstrips_z_sensor / ncopy;
325  offsetz = (ncopy % 2 == 0) ? -2. * cell_length_z / 2. * double(ncopy / 2) + cell_length_z / 2. + fphx_offset_z : -2. * cell_length_z / 2. * double(ncopy / 2) + fphx_offset_z;
326 
327  G4VPVParameterisation *fphxparam = new PHG4InttFPHXParameterisation(fphx_offsetx, +fphx_offsety, offsetz, 2. * cell_length_z / 2., ncopy);
328  new G4PVParameterised((boost::format("fphxcontainer_%d_%d") % inttlayer % itype).str(),
329  fphx_volume, fphxcontainer_volume, kZAxis, ncopy, fphxparam, OverlapCheck());
330 
331  // Glue for FPHX, silver powder epoxy, impletemented in the same way as FPHX
332  G4VSolid *fphx_glue_box = new G4Box((boost::format("fphx_glue_box_%d_%d") % inttlayer % itype).str(), fphx_glue_x / 2., fphx_y / 2., fphx_z / 2.);
333 
334  G4LogicalVolume *fphx_glue_volume = new G4LogicalVolume(fphx_glue_box, GetDetectorMaterial("SilverEpoxyGlue_INTT"),
335  (boost::format("fphx_glue_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
336  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
337  {
338  m_PassiveVolumeTuple.insert(std::make_pair(fphx_glue_volume, std::make_tuple(inttlayer, PHG4InttDefs::FPHX_GLUE)));
339  }
340  m_DisplayAction->AddVolume(fphx_glue_volume, "FPHXGlue");
341 
342  // Glue of FPHX Container
343  // make a container for the glue of FPHX chips, and then place them in the container
344  G4VSolid *fphx_gluecontainer_box = new G4Box((boost::format("fphx_gluecontainer_box_%d_%d") % inttlayer % itype).str(),
345  fphx_glue_x / 2., fphx_y / 2., hdi_z / 2.);
346  G4LogicalVolume *fphx_gluecontainer_volume = new G4LogicalVolume(fphx_gluecontainer_box, GetDetectorMaterial(rc->get_StringFlag("WorldMaterial")),
347  (boost::format("fphx_gluecontainer_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
348 
349  // Parameters for FPHX glue for G4VPVParameterisation are the same as FPGX's, so reuse them!
350  G4VPVParameterisation *fphx_glueparam = new PHG4InttFPHXParameterisation(fphx_offsetx, +fphx_offsety, offsetz, 2. * cell_length_z / 2., ncopy);
351 
352  new G4PVParameterised((boost::format("glue_fphxcontainer_%d_%d") % inttlayer % itype).str(),
353  fphx_glue_volume, fphx_gluecontainer_volume, kZAxis, ncopy, fphx_glueparam, OverlapCheck());
354  m_DisplayAction->AddVolume(fphx_gluecontainer_volume, "FPHXGlueContainer");
355 
356  double stave_x = 0.;
357  double stave_y = 0.;
358  G4LogicalVolume *stave_volume = NULL;
359  G4LogicalVolume *staveext_volume = NULL;
360 
361  // Carbon stave. This consists of the formed sheet, cooling water, the water tube, glue for the tube,
362  // rohacell foam to fill space around the tube, and the flat CFRP sheet, which completes the outer shell surrounds
363 
364  // Rohacel foam and cooling water pipe inside. Formed from straight sections and sections of a tube of
365  // radius 3.1905 mm. All have wall thickness of 0.1905 mm.
366  const double stave_thickness = params->get_double_param("stave_straight_cooler_x") * cm; // stave thickness
367  const double Rcmin = 0.30 * cm; // inner radius of curved section, same at both ends
368  const double Rcmax = Rcmin + stave_thickness; // outer radius of curved section, same at both ends
369  double Rcavge = (Rcmax + Rcmin) / 2.0;
370  double dphi_c = 23.19859051 * M_PI / 180.; // phi of the curved section
371  const double stave_z = hdi_z;
372 
373  // Make CFC structure
375  const double phic_begin[4] = {M_PI - dphi_c, -dphi_c, 0.0, M_PI};
376  const double dphic[4] = {dphi_c, dphi_c, dphi_c, dphi_c};
377 
378  G4Tubs *stave_curve_cons[4];
379  G4Tubs *stave_curve_ext_cons[4];
380  G4LogicalVolume *stave_curve_volume[4];
381  G4LogicalVolume *stave_curve_ext_volume[4];
382 
383  for (int i = 0; i < 4; i++)
384  {
385  stave_curve_cons[i] = new G4Tubs((boost::format("stave_curve_cons_%d_%d_%d") % inttlayer % itype % i).str(),
386  Rcmin, Rcmax, stave_z / 2., phic_begin[i], dphic[i]);
387  stave_curve_volume[i] = new G4LogicalVolume(stave_curve_cons[i], GetDetectorMaterial("CFRP_INTT"),
388  (boost::format("stave_curve_volume_%d_%d_%d") % inttlayer % itype % i).str(), 0, 0, 0);
389  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
390  {
391  m_PassiveVolumeTuple.insert(std::make_pair(stave_curve_volume[i], std::make_tuple(inttlayer, PHG4InttDefs::STAVE_CURVE)));
392  }
393  stave_curve_ext_cons[i] = new G4Tubs((boost::format("stave_curve_ext_cons_%d_%d_%d") % inttlayer % itype % i).str(),
394  Rcmin, Rcmax, hdiext_z / 2., phic_begin[i], dphic[i]);
395  stave_curve_ext_volume[i] = new G4LogicalVolume(stave_curve_ext_cons[i], GetDetectorMaterial("CFRP_INTT"),
396  (boost::format("stave_curve_ext_volume_%d_%d_%d") % inttlayer % itype % i).str(), 0, 0, 0);
397  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
398  {
399  m_PassiveVolumeTuple.insert(std::make_pair(stave_curve_ext_volume[i], std::make_tuple(inttlayer, PHG4InttDefs::STAVEEXT_CURVE)));
400  }
401  m_DisplayAction->AddVolume(stave_curve_volume[i], "StaveCurve");
402  m_DisplayAction->AddVolume(stave_curve_ext_volume[i], "StaveCurve");
403  }
404 
405  // we will need the length in y of the curved section as it is installed in the stave box
406  double curve_length_y = Rcavge * sin(dphi_c);
407 
408  // Make several straight sections for use in making the stave
409  double stave_straight_outer_y = params->get_double_param("stave_straight_outer_y") * cm;
410  double stave_straight_cooler_y = params->get_double_param("stave_straight_cooler_y") * cm;
411  double rohacell_straight_y = params->get_double_param("stave_straight_rohacell_y") * cm;
412 
413  // Outer straight sections of stave
414  G4VSolid *stave_straight_outer_box = new G4Box((boost::format("stave_straight_outer_box_%d_%d") % inttlayer % itype).str(),
415  stave_thickness / 2., stave_straight_outer_y / 2., stave_z / 2.);
416  G4LogicalVolume *stave_straight_outer_volume = new G4LogicalVolume(stave_straight_outer_box, GetDetectorMaterial("CFRP_INTT"),
417  (boost::format("stave_straight_outer_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
418  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
419  {
420  m_PassiveVolumeTuple.insert(std::make_pair(stave_straight_outer_volume, std::make_tuple(inttlayer, PHG4InttDefs::STAVE_STRAIGHT_OUTER)));
421  }
422  G4VSolid *stave_straight_outer_ext_box = new G4Box((boost::format("stave_straight_outer_ext_box_%d_%s") % inttlayer % itype).str(),
423  stave_thickness / 2., stave_straight_outer_y / 2., hdiext_z / 2.);
424  G4LogicalVolume *stave_straight_outer_ext_volume = new G4LogicalVolume(stave_straight_outer_ext_box, GetDetectorMaterial("CFRP_INTT"),
425  (boost::format("stave_straight_outer_ext_volume_%d_%s") % inttlayer % itype).str(), 0, 0, 0);
426  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
427  {
428  m_PassiveVolumeTuple.insert(std::make_pair(stave_straight_outer_ext_volume, std::make_tuple(inttlayer, PHG4InttDefs::STAVEEXT_STRAIGHT_OUTER)));
429  }
430 
431  //Top surface of stave
432  G4VSolid *stave_straight_cooler_box = new G4Box((boost::format("stave_straight_cooler_box_%d_%d") % inttlayer % itype).str(),
433  stave_thickness / 2., stave_straight_cooler_y / 2., stave_z / 2.);
434  G4LogicalVolume *stave_straight_cooler_volume = new G4LogicalVolume(stave_straight_cooler_box, GetDetectorMaterial("CFRP_INTT"),
435  (boost::format("stave_straight_cooler_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
436  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
437  {
438  m_PassiveVolumeTuple.insert(std::make_pair(stave_straight_cooler_volume, std::make_tuple(inttlayer, PHG4InttDefs::STAVE_STRAIGHT_COOLER)));
439  }
440  G4VSolid *stave_straight_cooler_ext_box = new G4Box((boost::format("stave_straight_cooler_ext_box_%d_%d") % inttlayer % itype).str(),
441  stave_thickness / 2., stave_straight_cooler_y / 2., hdiext_z / 2.);
442  G4LogicalVolume *stave_straight_cooler_ext_volume = new G4LogicalVolume(stave_straight_cooler_ext_box, GetDetectorMaterial("CFRP_INTT"),
443  (boost::format("stave_straight_cooler_ext_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
444  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
445  {
446  m_PassiveVolumeTuple.insert(std::make_pair(stave_straight_cooler_ext_volume, std::make_tuple(inttlayer, PHG4InttDefs::STAVEEXT_STRAIGHT_COOLER)));
447  }
448 
449  // Slant straight sections of stave
450  double stave_slant_cooler_y = params->get_double_param("stave_slant_cooler_y") * cm;
451  G4VSolid *stave_slant_cooler_box = new G4Box((boost::format("stave_slant_cooler_box_%d_%d") % inttlayer % itype).str(),
452  stave_thickness / 2., stave_slant_cooler_y / 2., stave_z / 2.);
453  G4LogicalVolume *stave_slant_cooler_volume = new G4LogicalVolume(stave_slant_cooler_box, GetDetectorMaterial("CFRP_INTT"),
454  (boost::format("stave_slant_cooler_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
455  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
456  {
457  m_PassiveVolumeTuple.insert(std::make_pair(stave_slant_cooler_volume, std::make_tuple(inttlayer, PHG4InttDefs::STAVE_STRAIGHT_COOLER)));
458  }
459  G4VSolid *stave_slant_cooler_ext_box = new G4Box((boost::format("stave_lant_cooler_ext_box_%d_%d") % inttlayer % itype).str(),
460  stave_thickness / 2., stave_slant_cooler_y / 2., hdiext_z / 2.);
461  G4LogicalVolume *stave_slant_cooler_ext_volume = new G4LogicalVolume(stave_slant_cooler_ext_box, GetDetectorMaterial("CFRP_INTT"),
462  (boost::format("stave_slant_cooler_ext_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
463  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
464  {
465  m_PassiveVolumeTuple.insert(std::make_pair(stave_slant_cooler_ext_volume, std::make_tuple(inttlayer, PHG4InttDefs::STAVEEXT_STRAIGHT_COOLER)));
466  }
467 
468  // Flat CFRP sheet on the bottom of the stave structure. It was introduced instead of PGS
469  G4VSolid *stave_bottom_cooler_box = new G4Box((boost::format("stave_bottom_cooler_box_%d_%d") % inttlayer % itype).str(),
470  stave_thickness / 2., hdi_y / 2., stave_z / 2.);
471 
472  G4LogicalVolume *stave_bottom_cooler_volume = new G4LogicalVolume(stave_bottom_cooler_box, GetDetectorMaterial("CFRP_INTT"),
473  (boost::format("stave_bottom_cooler_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
474  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
475  {
476  m_PassiveVolumeTuple.insert(std::make_pair(stave_bottom_cooler_volume, std::make_tuple(inttlayer, PHG4InttDefs::STAVE_BOTTOM_COOLER))); // should be changed soon
477  }
478 
479  G4VSolid *stave_bottom_cooler_ext_box = new G4Box((boost::format("stave_bottom_cooler_ext_box_%d_%s") % inttlayer % itype).str(), stave_thickness / 2., hdi_y / 2., hdiext_z / 2.);
480  G4LogicalVolume *stave_bottom_cooler_ext_volume = new G4LogicalVolume(stave_bottom_cooler_ext_box, GetDetectorMaterial("CFRP_INTT"),
481  (boost::format("stave_bottom_cooler_ext_volume_%d_%s") % inttlayer % itype).str(), 0, 0, 0);
482  if ((m_IsAbsorberActiveMap.find(inttlayer))->second > 0)
483  {
484  m_PassiveVolumeTuple.insert(std::make_pair(stave_bottom_cooler_ext_volume, std::make_tuple(inttlayer, PHG4InttDefs::STAVEEXT_BOTTOM_COOLER)));
485  }
486 
487  m_DisplayAction->AddVolume(stave_straight_cooler_volume, "StaveCooler");
488  m_DisplayAction->AddVolume(stave_straight_cooler_ext_volume, "StaveCooler");
489  m_DisplayAction->AddVolume(stave_straight_outer_volume, "StaveStraightOuter");
490  m_DisplayAction->AddVolume(stave_straight_outer_ext_volume, "StaveStraightOuter");
491  m_DisplayAction->AddVolume(stave_slant_cooler_volume, "StaveCooler");
492  m_DisplayAction->AddVolume(stave_slant_cooler_ext_volume, "StaveCooler");
493  m_DisplayAction->AddVolume(stave_bottom_cooler_volume, "StaveCooler");
494  m_DisplayAction->AddVolume(stave_bottom_cooler_ext_volume, "StaveCooler");
495 
496  // cooling pipe + water inside + glue outside
497  const double Rpmin = 0.10 * cm; // inner radius of cooling pipe section, same at both ends
498  const double Rpmax = 0.15 * cm; // outer radius of cooling pipe section, same at both ends
499  G4VSolid *stave_glue_box = new G4Box((boost::format("stave_glue_box_%d_%d") % inttlayer % itype).str(), 3. / 2, 3. / 2., stave_z / 2.);
500  G4LogicalVolume *stave_glue_volume = new G4LogicalVolume(stave_glue_box, GetDetectorMaterial("Epoxy"),
501  (boost::format("stave_glue_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
502  G4VSolid *staveext_glue_box = new G4Box((boost::format("staveext_glue_box_%d_%d") % inttlayer % itype).str(), 3. / 2., 3. / 2., hdiext_z / 2.);
503  G4LogicalVolume *staveext_glue_volume = new G4LogicalVolume(staveext_glue_box, GetDetectorMaterial("Epoxy"),
504  (boost::format("staveext_glue_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
505 
506  m_DisplayAction->AddVolume(stave_glue_volume, "StaveGlueBox");
507  m_DisplayAction->AddVolume(staveext_glue_volume, "StaveGlueBox");
508 
509  G4VSolid *stave_pipe_cons = new G4Tubs((boost::format("stave_pipe_cons_%d_%d") % inttlayer % itype).str(),
510  Rpmin, Rpmax, stave_z / 2., 0, 2.0 * M_PI);
511  G4LogicalVolume *stave_pipe_volume = new G4LogicalVolume(stave_pipe_cons, GetDetectorMaterial("CFRP_INTT"),
512  (boost::format("stave_pipe_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
513 
514  G4VSolid *staveext_pipe_cons = new G4Tubs((boost::format("staveext_pipe_cons_%d_%d") % inttlayer % itype).str(),
515  Rpmin, Rpmax, hdiext_z / 2., 0, 2.0 * M_PI);
516  G4LogicalVolume *staveext_pipe_volume = new G4LogicalVolume(staveext_pipe_cons, GetDetectorMaterial("CFRP_INTT"),
517  (boost::format("staveext_pipe_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
518 
519  m_DisplayAction->AddVolume(stave_pipe_volume, "StavePipe");
520  m_DisplayAction->AddVolume(staveext_pipe_volume, "StavePipe");
521 
522  G4VSolid *stave_water_cons = new G4Tubs((boost::format("stave_water_cons_%d_%d") % inttlayer % itype).str(),
523  0., Rpmin, stave_z / 2., 0, 2.0 * M_PI);
524  G4LogicalVolume *stave_water_volume = new G4LogicalVolume(stave_water_cons, GetDetectorMaterial("G4_WATER"),
525  (boost::format("stave_water_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
526 
527  G4VSolid *staveext_water_cons = new G4Tubs((boost::format("staveext_water_cons_%d_%d") % inttlayer % itype).str(),
528  0., Rpmin, hdiext_z / 2., 0, 2.0 * M_PI);
529  G4LogicalVolume *staveext_water_volume = new G4LogicalVolume(staveext_water_cons, GetDetectorMaterial("G4_WATER"),
530  (boost::format("staveext_water_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
531 
532  m_DisplayAction->AddVolume(stave_water_volume, "StaveWater");
533  m_DisplayAction->AddVolume(staveext_water_volume, "StaveWater");
534 
535  //rohacell foam
536  //straight boxes
537  G4VSolid *rohacell_straight_cons = new G4Box((boost::format("rohacell_straight_cons_%d_%d") % inttlayer % itype).str(), 3. / 2, rohacell_straight_y / 2., stave_z / 2.);
538  G4LogicalVolume *rohacell_straight_volume = new G4LogicalVolume(rohacell_straight_cons, GetDetectorMaterial("ROHACELL_FOAM_51"),
539  (boost::format("rohacell_straight_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
540 
541  G4VSolid *rohacellext_straight_cons = new G4Box((boost::format("rohacellext_straight_cons_%d_%d") % inttlayer % itype).str(), 3. / 2, rohacell_straight_y / 2., hdiext_z / 2.);
542  G4LogicalVolume *rohacellext_straight_volume = new G4LogicalVolume(rohacellext_straight_cons, GetDetectorMaterial("ROHACELL_FOAM_51"),
543  (boost::format("rohacellext_straight_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
544 
545  // make curved sections for rohacell foam
546  const double rh_phic_begin[2] = {-dphi_c, 0.0};
547  const double rh_dphic[2] = {dphi_c, dphi_c};
548  G4Tubs *rohacell_curve_cons[2];
549  G4LogicalVolume *rohacell_curve_volume[2];
550  G4Tubs *rohacellext_curve_cons[2];
551  G4LogicalVolume *rohacellext_curve_volume[2];
552  for (int i = 0; i < 2; i++)
553  {
554  rohacell_curve_cons[i] = new G4Tubs((boost::format("rohacell_curve_cons_%d_%d_%d") % inttlayer % itype % i).str(),
555  0., Rcmin, stave_z / 2., rh_phic_begin[i], rh_dphic[i]);
556  rohacell_curve_volume[i] = new G4LogicalVolume(rohacell_curve_cons[i], GetDetectorMaterial("ROHACELL_FOAM_51"),
557  (boost::format("rohacell_curve_volume_%d_%d_%d") % inttlayer % itype % i).str(), 0, 0, 0);
558  rohacellext_curve_cons[i] = new G4Tubs((boost::format("rohacellext_curve_cons_%d_%d_%d") % inttlayer % itype % i).str(),
559  0., Rcmin, hdiext_z / 2., rh_phic_begin[i], rh_dphic[i]);
560  rohacellext_curve_volume[i] = new G4LogicalVolume(rohacellext_curve_cons[i], GetDetectorMaterial("ROHACELL_FOAM_51"),
561  (boost::format("rohacellext_curve_volume_%d_%d_%d") % inttlayer % itype % i).str(), 0, 0, 0);
562  }
563 
564  // make trapezoidal sections for rohacell foam
565  G4GenericTrap *rohacell_trap_cons[2];
566  G4LogicalVolume *rohacell_trap_volume[2];
567  G4GenericTrap *rohacellext_trap_cons[2];
568  G4LogicalVolume *rohacellext_trap_volume[2];
569  for (int i = 0; i < 2; i++)
570  {
571  double shift = 1.e-5; // To mitigate fm order level overlaps reported by GEANT4...
572  std::vector<G4TwoVector> rohatrap(8);
573  if (i == 0)
574  {
575  rohatrap[0] = G4TwoVector(0. * cm, 0. * cm);
576  rohatrap[1] = G4TwoVector(Rcmin * cos(dphi_c) - shift, -Rcmin * sin(dphi_c));
577  rohatrap[2] = G4TwoVector(Rcmin * (1. - cos(dphi_c)) - shift, -stave_slant_cooler_y * cos(dphi_c) - Rcmin * sin(dphi_c));
578  rohatrap[3] = G4TwoVector(0. * cm, -stave_slant_cooler_y * cos(dphi_c) - Rcmin * sin(dphi_c));
579  }
580  else
581  {
582  rohatrap[0] = G4TwoVector(0. * cm, +stave_slant_cooler_y * cos(dphi_c) + Rcmin * sin(dphi_c));
583  rohatrap[1] = G4TwoVector(Rcmax * (1. - cos(dphi_c)) - shift, +stave_slant_cooler_y * cos(dphi_c) + Rcmin * sin(dphi_c));
584  rohatrap[2] = G4TwoVector(Rcmin * cos(dphi_c) - shift, +Rcmin * sin(dphi_c));
585  rohatrap[3] = G4TwoVector(0. * cm, 0. * cm);
586  }
587  rohatrap[4] = rohatrap[0];
588  rohatrap[5] = rohatrap[1];
589  rohatrap[6] = rohatrap[2];
590  rohatrap[7] = rohatrap[3];
591 
592  rohacell_trap_cons[i] = new G4GenericTrap((boost::format("rohacell_trap_cons_%d_%d_%d") % inttlayer % itype % i).str(), stave_z / 2., rohatrap);
593  rohacell_trap_volume[i] = new G4LogicalVolume(rohacell_trap_cons[i], GetDetectorMaterial("ROHACELL_FOAM_51"),
594  (boost::format("rohacell_trap_volume_%d_%d_%d") % inttlayer % itype % i).str(), 0, 0, 0);
595 
596  rohacellext_trap_cons[i] = new G4GenericTrap((boost::format("rohacellext_trap_cons_%d_%d_%d") % inttlayer % itype % i).str(), hdiext_z / 2., rohatrap);
597  rohacellext_trap_volume[i] = new G4LogicalVolume(rohacellext_trap_cons[i], GetDetectorMaterial("ROHACELL_FOAM_51"),
598  (boost::format("rohacellext_trap_volume_%d_%d_%d") % inttlayer % itype % i).str(), 0, 0, 0);
599  }
600 
601  m_DisplayAction->AddVolume(rohacell_straight_volume, "RohaCell");
602  m_DisplayAction->AddVolume(rohacellext_straight_volume, "RohaCell");
603  for (int i = 0; i < 2; i++)
604  {
605  m_DisplayAction->AddVolume(rohacell_curve_volume[i], "RohaCell");
606  m_DisplayAction->AddVolume(rohacellext_curve_volume[i], "RohaCell");
607  m_DisplayAction->AddVolume(rohacell_trap_volume[i], "RohaCell");
608  m_DisplayAction->AddVolume(rohacellext_trap_volume[i], "RohaCell");
609  }
610 
611  // Now we combine the elements of a stave defined above into a stave
612  // Create a stave volume to install the stave sections into. The volume has to be big enouigh to contain the cooling tube
613  double cooler_gap_x = 0.3 * cm; // id of cooling tube in cm
614  double cooler_wall = stave_thickness; // outer wall thickness of cooling tube
615  double cooler_x = cooler_gap_x + 2.0 * cooler_wall; // thickness of the formed sheet, the flat sheet, and the gap b/w the sheets
616  stave_x = cooler_x;
617  stave_y = hdi_y;
618 
619  // Make stave volume. Drop two corners in positive x to prevent ladder_volume overlapping
620  // with neighbouring ladders because of small clearance in the latest configuration
621  G4RotationMatrix *stv_rot_pos = new G4RotationMatrix();
622  stv_rot_pos->rotateZ(-15. * M_PI / 180.);
623  G4ThreeVector stvTranspos(stave_x / 2., stave_y / 2., 0.);
624 
625  G4RotationMatrix *stv_rot_neg = new G4RotationMatrix();
626  stv_rot_neg->rotateZ(+15. * M_PI / 180.);
627  G4ThreeVector stvTransneg(stave_x / 2., -stave_y / 2., 0.);
628 
629  G4VSolid *stave_basebox = new G4Box((boost::format("stave_basebox_%d_%d") % inttlayer % itype).str(), stave_x / 2., stave_y / 2., stave_z / 2.);
630  G4VSolid *stave_subtbox = new G4Box((boost::format("stave_subtbox_%d_%d") % inttlayer % itype).str(), stave_x / 1.5, stave_y / 1.5, stave_z / 1.); // has to be longer in z to avoid coincident surface
631 
632  G4VSolid *stave_box1 = new G4SubtractionSolid((boost::format("stave_box1_%d_%d") % inttlayer % itype).str(), stave_basebox, stave_subtbox, stv_rot_pos, stvTranspos);
633 
634  G4VSolid *stave_box = new G4SubtractionSolid((boost::format("stave_box_%d_%d") % inttlayer % itype).str(), stave_box1, stave_subtbox, stv_rot_neg, stvTransneg);
635 
636  stave_volume = new G4LogicalVolume(stave_box, GetDetectorMaterial(rc->get_StringFlag("WorldMaterial")),
637  (boost::format("stave_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
638 
639  G4VSolid *staveext_basebox = new G4Box((boost::format("staveext_basebox_%d_%d") % inttlayer % itype).str(), stave_x / 2., stave_y / 2., hdiext_z / 2.);
640  G4VSolid *staveext_subtbox = new G4Box((boost::format("staveext_subtbox_%d_%d") % inttlayer % itype).str(), stave_x / 1.5, stave_y / 1.5, hdiext_z / 1.); // has to be longer in z to avoid coincident surface
641 
642  G4VSolid *staveext_box1 = new G4SubtractionSolid((boost::format("staveext_box1_%d_%d") % inttlayer % itype).str(), staveext_basebox, staveext_subtbox, stv_rot_pos, stvTranspos);
643 
644  G4VSolid *staveext_box = new G4SubtractionSolid((boost::format("staveext_box_%d_%d") % inttlayer % itype).str(), staveext_box1, staveext_subtbox, stv_rot_neg, stvTransneg);
645 
646  staveext_volume = new G4LogicalVolume(staveext_box, GetDetectorMaterial(rc->get_StringFlag("WorldMaterial")),
647  (boost::format("staveext_volume_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
648  // the rotation matrices are just used by G4VSolid, ownership is not taken over
649  delete stv_rot_pos;
650  delete stv_rot_neg;
651 
652  m_DisplayAction->AddVolume(stave_volume, "StaveBox");
653  m_DisplayAction->AddVolume(staveext_volume, "StaveBox");
654 
655  // Assemble the elements into the stave volume and the stave extension volume
656  // They are place relative to the center of the stave box. Thus the offset of the center of the segment is relative to the center of the satev box.
657  // But we want the segment to be located relative to the lowest x limit of the stave box.
658  if (laddertype == PHG4InttDefs::SEGMENTATION_Z) // Obsolete!!
659  {
660  // only one cooling tube in laddertype 0
661  // Place the straight sections. We add the middle, then above x axis, then below x axis
662  double x_off_str[3] =
663  {
664  Rcavge - stave_x / 2.,
665  (Rcmax - Rcmin) / 2. - stave_x / 2.,
666  (Rcmax - Rcmin) / 2. - stave_x / 2.};
667  double y_off_str[3] =
668  {
669  0.0,
670  +stave_straight_cooler_y / 2. + 2. * curve_length_y + stave_straight_outer_y / 2.,
671  -stave_straight_cooler_y / 2. - 2. * curve_length_y - stave_straight_outer_y / 2.};
672 
673  for (int i = 0; i < 3; i++)
674  {
675  if (i == 0)
676  {
677  new G4PVPlacement(0, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0), stave_straight_cooler_volume,
678  (boost::format("stave_straight_cooler_%d_%d_%d") % i % inttlayer % itype).str(), stave_volume, false, 0, OverlapCheck());
679  new G4PVPlacement(0, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0), stave_straight_cooler_ext_volume,
680  (boost::format("stave_straight_cooler_ext_%d_%d_%d") % i % inttlayer % itype).str(), staveext_volume, false, 0, OverlapCheck());
681  }
682  else
683  {
684  new G4PVPlacement(0, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0), stave_straight_outer_volume,
685  (boost::format("stave_straight_outer_%d_%d_%d") % i % inttlayer % itype).str(), stave_volume, false, 0, OverlapCheck());
686  new G4PVPlacement(0, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0), stave_straight_outer_ext_volume,
687  (boost::format("stave_straight_outer_ext_%d_%d_%d") % i % inttlayer % itype).str(), staveext_volume, false, 0, OverlapCheck());
688  }
689  }
690  // The cooler curved sections are made using 2 curved sections in a recurve on each side of the cooler straight section
691  // The tube sections used here have the origin of their volume at their center of rotation. Rcavge
692  // Each curve section is moved to the center of the stave volume by a translation of +/- Rcavge
693  // Then it is moved to the outside or the inside of the stave volume by a translation of +/- cooler_gap_x / 2.
694  // we start at lowest y and work up in y
695 
696  double x_off_cooler[4] =
697  {
698  Rcavge - cooler_gap_x / 2.,
699  -Rcavge + cooler_gap_x / 2.,
700  -Rcavge + cooler_gap_x / 2.,
701  Rcavge - cooler_gap_x / 2.};
702  double y_off_cooler[4] =
703  {
704  -stave_straight_cooler_y / 2. - 2. * curve_length_y,
705  -stave_straight_cooler_y / 2.,
706  +stave_straight_cooler_y / 2.,
707  +stave_straight_cooler_y / 2. + 2. * curve_length_y};
708 
709  for (int i = 0; i < 4; i++)
710  {
711  new G4PVPlacement(0, G4ThreeVector(x_off_cooler[i], y_off_cooler[i], 0.0), stave_curve_volume[i],
712  (boost::format("stave_curve_%d_%d_%d") % inttlayer % itype % i).str(), stave_volume, false, 0, OverlapCheck());
713  new G4PVPlacement(0, G4ThreeVector(x_off_cooler[i], y_off_cooler[i], 0.0), stave_curve_ext_volume[i],
714  (boost::format("stave_curve_ext_%d_%d_%s") % inttlayer % itype % i).str(), staveext_volume, false, 0, OverlapCheck());
715  }
716  }
717  else if (laddertype == PHG4InttDefs::SEGMENTATION_PHI) // The type PHG4InttDefs::SEGMENTATION_PHI ladder
718  {
719  // First place the straight sections, do the extension at the same time
720  // we alternate positive and negative y values here
721  double x_off_str[6] =
722  {
723  (Rcmax + Rcmin) / 2. - stave_x / 2. + stave_thickness, // inner straight section
724  (Rcmax + Rcmax) / 4. - stave_x / 2. + stave_thickness, // slant section
725  (Rcmax + Rcmax) / 4. - stave_x / 2. + stave_thickness, // slant section
726  (Rcmax - Rcmin) / 2. - stave_x / 2. + stave_thickness, // outer straight section
727  (Rcmax - Rcmin) / 2. - stave_x / 2. + stave_thickness, // outer straight section
728  (Rcmax - Rcmin) / 2. - stave_x / 2. // bottom section
729  };
730  double y_off_str[6] =
731  {
732  0.0, // inner straight section
733  -stave_straight_cooler_y / 2. - 1. * curve_length_y - cos(dphi_c) * stave_slant_cooler_y / 2., // slant section
734  +stave_straight_cooler_y / 2. + 1. * curve_length_y + cos(dphi_c) * stave_slant_cooler_y / 2., // slant section
735  -stave_straight_cooler_y / 2. - 2. * curve_length_y - cos(dphi_c) * stave_slant_cooler_y - stave_straight_outer_y / 2., // outer straight section
736  +stave_straight_cooler_y / 2. + 2. * curve_length_y + cos(dphi_c) * stave_slant_cooler_y + stave_straight_outer_y / 2., // outer straight section
737  0.0
738  // bottom straight section
739  };
740 
741  for (int i = 0; i < 6; i++)
742  {
743  if (i == 0) // inner straight section
744  {
745  new G4PVPlacement(0, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0), stave_straight_cooler_volume,
746  (boost::format("stave_straight_cooler_%d_%d_%d") % inttlayer % itype % i).str(), stave_volume, false, 0, OverlapCheck());
747  new G4PVPlacement(0, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0), stave_straight_cooler_ext_volume,
748  (boost::format("stave_straight_cooler_ext_%d_%d_%s") % inttlayer % itype % i).str(), staveext_volume, false, 0, OverlapCheck());
749  }
750  else if (i == 1 || i == 2) // slant section
751  {
752  G4RotationMatrix rotation;
753  if (i == 1)
754  rotation.rotateZ(-1. * dphi_c);
755  else if (i == 2)
756  rotation.rotateZ(dphi_c);
757  new G4PVPlacement(G4Transform3D(rotation, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0)), stave_slant_cooler_volume,
758  (boost::format("stave_slant_cooler_%d_%d_%d") % inttlayer % itype % i).str(), stave_volume, false, 0, OverlapCheck());
759  new G4PVPlacement(G4Transform3D(rotation, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0)), stave_slant_cooler_ext_volume,
760  (boost::format("stave_slant_cooler_ext_%d_%d_%d") % inttlayer % itype % i).str(), staveext_volume, false, 0, OverlapCheck());
761  }
762  else if (i == 3 || i == 4) // outer straight section
763  {
764  new G4PVPlacement(0, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0), stave_straight_outer_volume,
765  (boost::format("stave_straight_outer_%d_%d_%d") % inttlayer % itype % i).str(), stave_volume, false, 0, OverlapCheck());
766  new G4PVPlacement(0, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0), stave_straight_outer_ext_volume,
767  (boost::format("stave_straight_outer_ext_%d_%d_%s") % inttlayer % itype % i).str(), staveext_volume, false, 0, OverlapCheck());
768  }
769  else // bottom
770  {
771  new G4PVPlacement(0, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0), stave_bottom_cooler_volume,
772  (boost::format("stave_bottom_cooler_%d_%d_%d") % inttlayer % itype % i).str(), stave_volume, false, 0, OverlapCheck());
773  new G4PVPlacement(0, G4ThreeVector(x_off_str[i], y_off_str[i], 0.0), stave_bottom_cooler_ext_volume,
774  (boost::format("stave_bottom_cooler_ext_%d_%d_%s") % inttlayer % itype % i).str(), staveext_volume, false, 0, OverlapCheck());
775  }
776  }
777 
780 
781  double x_off_curve[4] =
782  {
783  // increasing in y
784  +Rcavge - cooler_gap_x / 2. + stave_thickness / 2.,
785  -Rcavge + cooler_gap_x / 2. + stave_thickness / 2.,
786  -Rcavge + cooler_gap_x / 2. + stave_thickness / 2.,
787  +Rcavge - cooler_gap_x / 2. + stave_thickness / 2.};
788  double y_off_curve[4] =
789  {
790  // increasing in y
791  -stave_straight_cooler_y / 2. - 2. * curve_length_y - cos(dphi_c) * stave_slant_cooler_y,
792  -stave_straight_cooler_y / 2.,
793  +stave_straight_cooler_y / 2.,
794  +stave_straight_cooler_y / 2. + 2. * curve_length_y + cos(dphi_c) * stave_slant_cooler_y};
795 
796  for (int i = 0; i < 4; i++)
797  {
798  new G4PVPlacement(0, G4ThreeVector(x_off_curve[i], y_off_curve[i], 0.0), stave_curve_volume[i], (boost::format("stave_curve_%d_%d_%d") % inttlayer % itype % i).str(), stave_volume, false, 0, OverlapCheck());
799  new G4PVPlacement(0, G4ThreeVector(x_off_curve[i], y_off_curve[i], 0.0), stave_curve_ext_volume[i], (boost::format("stave_curve_ext_%d_%d_%s") % inttlayer % itype % i).str(), staveext_volume, false, 0, OverlapCheck());
800  }
801 
802  // Place the rohacell foam
803  // straight box section
804  double x_off_roha_str[2] =
805  {
806  // increasing in y
807  -cooler_wall / 2. + stave_thickness / 2.,
808  -cooler_wall / 2. + stave_thickness / 2.};
809  double y_off_roha_str[2] =
810  {
811  // increasing in y
812  -3. / 2. - rohacell_straight_y / 2.,
813  +3. / 2. + rohacell_straight_y / 2.};
814 
815  for (int i = 0; i < 2; i++)
816  {
817  new G4PVPlacement(0, G4ThreeVector(x_off_roha_str[i], y_off_roha_str[i], 0.0), rohacell_straight_volume, (boost::format("rohacell_straight_%d_%d_%d") % inttlayer % itype % i).str(), stave_volume, false, 0, OverlapCheck());
818  new G4PVPlacement(0, G4ThreeVector(x_off_roha_str[i], y_off_roha_str[i], 0.0), rohacellext_straight_volume, (boost::format("rohacell_straight_ext_%d_%d_%d") % inttlayer % itype % i).str(), staveext_volume, false, 0, OverlapCheck());
819  }
820 
822  double x_off_roha_curve[2] =
823  {
824  // increasing in y
825  -Rcavge + cooler_gap_x / 2. + stave_thickness / 2.,
826  -Rcavge + cooler_gap_x / 2. + stave_thickness / 2.};
827  double y_off_roha_curve[2] =
828  {
829  // increasing in y
830  -3. / 2. - rohacell_straight_y,
831  +3. / 2. + rohacell_straight_y};
832 
833  for (int i = 0; i < 2; i++)
834  {
835  new G4PVPlacement(0, G4ThreeVector(x_off_roha_curve[i], y_off_roha_curve[i], 0.0), rohacell_curve_volume[i], (boost::format("rohacell_curve_%d_%d_%d") % inttlayer % itype % i).str(), stave_volume, false, 0, OverlapCheck());
836  new G4PVPlacement(0, G4ThreeVector(x_off_roha_curve[i], y_off_roha_curve[i], 0.0), rohacellext_curve_volume[i], (boost::format("rohacell_curve_ext_%d_%d_%d") % inttlayer % itype % i).str(), staveext_volume, false, 0, OverlapCheck());
837  }
838 
839  // trapezoidal section
840  double x_off_roha_trap[2] =
841  {
842  // increasing in y
843  -Rcmin - cooler_wall / 2. + cooler_gap_x / 2. + stave_thickness / 2.,
844  -Rcmin - cooler_wall / 2. + cooler_gap_x / 2. + stave_thickness / 2.};
845  double y_off_roha_trap[2] =
846  {
847  // increasing in y
848  -3. / 2. - rohacell_straight_y,
849  +3. / 2. + rohacell_straight_y};
850 
851  for (int i = 0; i < 2; i++)
852  {
853  new G4PVPlacement(0, G4ThreeVector(x_off_roha_trap[i], y_off_roha_trap[i], 0.0), rohacell_trap_volume[i], (boost::format("rohacell_trap_%d_%d_%d") % inttlayer % itype % i).str(), stave_volume, false, 0, OverlapCheck());
854  new G4PVPlacement(0, G4ThreeVector(x_off_roha_trap[i], y_off_roha_trap[i], 0.0), rohacellext_trap_volume[i], (boost::format("rohacell_trap_ext_%d_%d_%d") % inttlayer % itype % i).str(), staveext_volume, false, 0, OverlapCheck());
855  }
856 
857  // place glue box, cooling pipe and water inside
858  new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, 0.0), stave_pipe_volume, (boost::format("stave_pipe_%d_%d") % inttlayer % itype).str(), stave_glue_volume, false, 0, OverlapCheck());
859  new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, 0.0), staveext_pipe_volume, (boost::format("stave_pipe_ext_%d_%d") % inttlayer % itype).str(), staveext_glue_volume, false, 0, OverlapCheck());
860  new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, 0.0), stave_water_volume, (boost::format("stave_water_%d_%d") % inttlayer % itype).str(), stave_glue_volume, false, 0, OverlapCheck());
861  new G4PVPlacement(0, G4ThreeVector(0.0, 0.0, 0.0), staveext_water_volume, (boost::format("stave_water_ext_%d_%d") % inttlayer % itype).str(), staveext_glue_volume, false, 0, OverlapCheck());
862 
863  // place of stave_glue_volume -cooler_wall / 2. + stave_thickness / 2. is actually 0. But I don't put 0 directry to make the origin of the value clear
864  new G4PVPlacement(0, G4ThreeVector(-cooler_wall / 2. + stave_thickness / 2., 0.0, 0.0), stave_glue_volume, (boost::format("stave_glue_%d_%d") % inttlayer % itype).str(), stave_volume, false, 0, OverlapCheck());
865  new G4PVPlacement(0, G4ThreeVector(-cooler_wall / 2. + stave_thickness / 2., 0.0, 0.0), staveext_glue_volume, (boost::format("stave_glue_ext_%d_%d") % inttlayer % itype).str(), staveext_volume, false, 0, OverlapCheck());
866  }
867  else
868  {
869  std::cout << PHWHERE << "invalid laddertype " << laddertype << std::endl;
870  gSystem->Exit(1);
871  }
872 
873  // ----- Step 2 ---------------------------------------------------
874  // We place Si-sensor, FPHX, HDI, and stave in the ladder volume.
875  // ================================================================
876 
877  // Make the ladder volume first
878  // We are still in the loop over inner or outer sensors. This is the ladder volume corresponding to this sensor.
879  // But the thickness of the glue for FPHX is used since it's taller than the sensor in x.
880  const double ladder_x = stave_x + hdi_kapton_x + hdi_copper_x + fphx_glue_x + fphx_x;
881  double ladder_y = hdi_y;
882 
883  // Make ladder volume. Drop two corners in positive x as done for stave volume
884  G4RotationMatrix *lad_box_rotpos = new G4RotationMatrix();
885  lad_box_rotpos->rotateZ(-15. * M_PI / 180.);
886  G4ThreeVector ladTranspos(ladder_x / 2., ladder_y / 2., 0.);
887 
888  G4RotationMatrix *lad_box_rotneg = new G4RotationMatrix();
889  lad_box_rotneg->rotateZ(+15. * M_PI / 180.);
890  G4ThreeVector ladTransneg(ladder_x / 2., -ladder_y / 2., 0.);
891 
892  G4VSolid *ladder_basebox = new G4Box((boost::format("ladder_basebox_%d_%d") % inttlayer % itype).str(), ladder_x / 2., ladder_y / 2., hdi_z / 2.);
893  G4VSolid *ladder_subtbox = new G4Box((boost::format("ladder_subtbox_%d_%d") % inttlayer % itype).str(), stave_x / 1.5, ladder_y / 1.5, hdi_z / 1.); // has to be longer in z to avoid coincident surface
894 
895  G4VSolid *ladder_box1 = new G4SubtractionSolid((boost::format("ladder_box1_%d_%d") % inttlayer % itype).str(), ladder_basebox, ladder_subtbox, lad_box_rotpos, ladTranspos);
896 
897  G4VSolid *ladder_box = new G4SubtractionSolid((boost::format("ladder_box_%d_%d") % inttlayer % itype).str(), ladder_box1, ladder_subtbox, lad_box_rotneg, ladTransneg);
898 
899  G4LogicalVolume *ladder_volume = new G4LogicalVolume(ladder_box, GetDetectorMaterial(rc->get_StringFlag("WorldMaterial")), (boost::format("ladder_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
900 
901  G4VSolid *ladderext_basebox = new G4Box((boost::format("ladderext_basebox_%d_%d") % inttlayer % itype).str(), ladder_x / 2., ladder_y / 2., hdiext_z / 2.);
902  G4VSolid *ladderext_subtbox = new G4Box((boost::format("ladderext_subtbox_%d_%d") % inttlayer % itype).str(), stave_x / 1.5, ladder_y / 1.5, hdiext_z / 1.); // has to be longer in z to avoid coincident surface
903 
904  G4VSolid *ladderext_box1 = new G4SubtractionSolid((boost::format("ladderext_box1_%d_%d") % inttlayer % itype).str(), ladderext_basebox, ladderext_subtbox, lad_box_rotpos, ladTranspos);
905 
906  G4VSolid *ladderext_box = new G4SubtractionSolid((boost::format("ladderext_box_%d_%d") % inttlayer % itype).str(), ladderext_box1, ladderext_subtbox, lad_box_rotneg, ladTransneg);
907 
908  G4LogicalVolume *ladderext_volume = new G4LogicalVolume(ladderext_box, GetDetectorMaterial(rc->get_StringFlag("WorldMaterial")), (boost::format("ladderext_%d_%d") % inttlayer % itype).str(), 0, 0, 0);
909  // the rotation matrices are just used by G4VSolid, ownership is not taken over
910  delete lad_box_rotpos;
911  delete lad_box_rotneg;
912  m_DisplayAction->AddVolume(ladder_volume, "Ladder");
913  m_DisplayAction->AddVolume(ladderext_volume, "Ladder");
914 
915  // Now add the components of the ladder to the ladder volume
916  // The sensor is closest to the beam pipe, the stave cooler is furthest away
917  // Note that the cooler has been assembled in the stave volume with the top at larger x, so the sensor will be at smaller x
918  // That will be the configuration when the ladder is at phi = 0 degrees, the positive x direction
919 
920  // We start at the most positive x value and add the stave first
921 
922  // Carbon stave
923  double TVstave_y = 0.0;
924  const double TVstave_x = ladder_x / 2. - stave_x / 2.;
925  new G4PVPlacement(0, G4ThreeVector(TVstave_x, TVstave_y, 0.0), stave_volume, (boost::format("stave_%d_%d") % inttlayer % itype).str(),
926  ladder_volume, false, 0, OverlapCheck());
927  new G4PVPlacement(0, G4ThreeVector(TVstave_x, TVstave_y, 0.0), staveext_volume, (boost::format("staveext_%d_%s") % inttlayer % itype).str(),
928  ladderext_volume, false, 0, OverlapCheck());
929 
930  // HDI Kapton
931  const double TVhdi_kapton_x = TVstave_x - stave_x / 2. - hdi_kapton_x / 2.;
932  new G4PVPlacement(0, G4ThreeVector(TVhdi_kapton_x, TVstave_y, 0.0), hdi_kapton_volume, (boost::format("hdikapton_%d_%d") % inttlayer % itype).str(), ladder_volume, false, 0, OverlapCheck());
933  new G4PVPlacement(0, G4ThreeVector(TVhdi_kapton_x, TVstave_y, 0.0), hdiext_kapton_volume, (boost::format("hdiextkapton_%d_%s") % inttlayer % itype).str(), ladderext_volume, false, 0, OverlapCheck());
934 
935  // HDI copper
936  const double TVhdi_copper_x = TVhdi_kapton_x - hdi_kapton_x / 2. - hdi_copper_x / 2.;
937  new G4PVPlacement(0, G4ThreeVector(TVhdi_copper_x, TVstave_y, 0.0), hdi_copper_volume, (boost::format("hdicopper_%d_%d") % inttlayer % itype).str(), ladder_volume, false, 0, OverlapCheck());
938  new G4PVPlacement(0, G4ThreeVector(TVhdi_copper_x, TVstave_y, 0.0), hdiext_copper_volume, (boost::format("hdiextcopper_%d_%s") % inttlayer % itype).str(), ladderext_volume, false, 0, OverlapCheck());
939 
940  // Glue for Si-sensor
941  const double TVsi_glue_x = TVhdi_copper_x - hdi_copper_x / 2. - si_glue_x / 2.;
942  new G4PVPlacement(0, G4ThreeVector(TVsi_glue_x, TVstave_y, 0.0), si_glue_volume, (boost::format("si_glue_%d_%d") % inttlayer % itype).str(), ladder_volume, false, 0, OverlapCheck());
943 
944  // Si-sensor
945  double TVSi_y = 0.0;
946  // sensor is not centered in y in the ladder volume for the Z sensitive ladders
947  if (laddertype == PHG4InttDefs::SEGMENTATION_Z)
948  TVSi_y = +sensor_offset_y;
949  //const double TVSi_x = TVhdi_copper_x - hdi_copper_x / 2. - siactive_x / 2.;
950  const double TVSi_x = TVsi_glue_x - si_glue_x / 2. - siactive_x / 2.;
951  new G4PVPlacement(0, G4ThreeVector(TVSi_x, TVSi_y, 0.0), siinactive_volume,
952  (boost::format("siinactive_%d_%d") % inttlayer % itype).str(), ladder_volume, false, 0, OverlapCheck());
953  new G4PVPlacement(0, G4ThreeVector(TVSi_x, TVSi_y, 0.0), siactive_volume,
954  (boost::format("siactive_%d_%d") % inttlayer % itype).str(), ladder_volume, false, 0, OverlapCheck());
955 
956  // FPHX glue
957  const double TVfphx_glue_x = TVhdi_copper_x - hdi_copper_x / 2. - fphx_glue_x / 2.;
958  double TVfphx_glue_y = sifull_y / 2. + gap_sensor_fphx + fphx_y / 2.;
959  if (laddertype == PHG4InttDefs::SEGMENTATION_Z)
960  TVfphx_glue_y -= sensor_offset_y;
961 
962  // laddertype PHG4InttDefs::SEGMENTATION_Z has only one FPHX, laddertype PHG4InttDefs::SEGMENTATION_PHI has two
963  if (laddertype == PHG4InttDefs::SEGMENTATION_PHI)
964  {
965  new G4PVPlacement(0, G4ThreeVector(TVfphx_glue_x, +TVfphx_glue_y, 0.0), fphx_gluecontainer_volume, (boost::format("fphx_gluecontainerp_%d_%d") % inttlayer % itype).str(), ladder_volume, false, 0, OverlapCheck());
966  }
967  new G4PVPlacement(0, G4ThreeVector(TVfphx_glue_x, -TVfphx_glue_y, 0.0), fphx_gluecontainer_volume, (boost::format("fphx_gluecontainerm_%d_%d") % inttlayer % itype).str(), ladder_volume, false, 0, OverlapCheck());
968 
969  // FPHX container
970  const double TVfphx_x = TVfphx_glue_x - fphx_glue_x / 2. - fphx_x / 2.;
971  double TVfphx_y = sifull_y / 2. + gap_sensor_fphx + fphx_y / 2.;
972  if (laddertype == PHG4InttDefs::SEGMENTATION_Z)
973  TVfphx_y -= sensor_offset_y;
974 
975  // laddertype PHG4InttDefs::SEGMENTATION_Z has only one FPHX, laddertype PHG4InttDefs::SEGMENTATION_PHI has two
976  if (laddertype == PHG4InttDefs::SEGMENTATION_PHI)
977  {
978  new G4PVPlacement(0, G4ThreeVector(TVfphx_x, +TVfphx_y, 0.0), fphxcontainer_volume, (boost::format("fphxcontainerp_%d_%d") % inttlayer % itype).str(), ladder_volume, false, 0, OverlapCheck());
979  }
980  new G4PVPlacement(0, G4ThreeVector(TVfphx_x, -TVfphx_y, 0.0), fphxcontainer_volume, (boost::format("fphxcontainerm_%d_%d") % inttlayer % itype).str(), ladder_volume, false, 0, OverlapCheck());
981 
982  // ----- Step 3 -----------------------------------------------------------------------------------------------
983  // We install the section of ladder for this sensor at all requested phi values and at positive and negative Z
984  //=============================================================================================================
985 
986  // Distribute Ladders in phi
987  // We are still in the loops over layer and sensor type, we will place copies of the ladder section for this sensor
988  // at all ladder phi values, and at both positive and negative Z values.
989 
990  // given radius values are for the center of the sensor, we need the x offset from center of ladder to center of sensor so we can place the ladder
991  double sensor_offset_x_ladder = 0.0 - TVSi_x; // ladder center is at x = 0.0 by construction. Sensor is at lower x, so TVSi_x is negative
992 
993  const double dphi = 2 * M_PI / nladders_layer;
994 
995  m_PosZ[inttlayer][itype] = (itype == 0) ? hdi_z / 2. : hdi_z_arr[inttlayer][0] + hdi_z / 2.; // location of center of ladder in Z
996  m_StripOffsetX[inttlayer] = sensor_offset_x_ladder;
997 
998  // The sensors have no tilt in the new design
999  // The type 1 ladders have the sensor at the center of the ladder in phi, so that is easy
1000  // The type 0 ladders are more complicated because the sensor center is perpendicular to the radial vector and the sensor is not at the ladder center
1001  // We made the stave box symmetric in y around the sensor center to simplify things
1002 
1003  for (int icopy = 0; icopy < nladders_layer; icopy++)
1004  {
1005  // sensor center
1006  const double phi = offsetphi + dphi * icopy; // if offsetphi is zero we start at zero
1007 
1008  double radius;
1009  // Make each layer at a single radius - i.e. what was formerly a sub-layer is now considered a layer
1010  radius = m_SensorRadius[inttlayer];
1011  radius += sensor_offset_x_ladder;
1012 
1013  double p = 0.0;
1014  if (laddertype == PHG4InttDefs::SEGMENTATION_Z)
1015  {
1016  // The Z sensitive ladders have the sensors offset in y relative to the ladder center
1017  // We have to slightly rotate the ladder in its own frame to make the radial vector to the sensor center normal to the sensor face
1018  p = atan(sensor_offset_y / radius);
1019  // then we adjust the distance to the center of the ladder to put the sensor at the requested distance from the center of the barrel
1020  radius /= cos(p);
1021  }
1022 
1023  // these describe the center of the ladder volume, placing it so that the center of the sensor is at phi = dphi * icopy, and at the correct radius
1024  const double posx = radius * cos(phi - p);
1025  const double posy = radius * sin(phi - p);
1026  const double fRotate = p + (phi - p) + offsetrot; // rotate in its own frame to make sensor perp to radial vector (p), then additionally rotate to account for ladder phi
1027  G4RotationMatrix ladderrotation;
1028  ladderrotation.rotateZ(fRotate);
1029 
1030  // this placement version rotates the ladder in its own frame by fRotate, then translates the center to posx, posy, +/- m_PosZ
1031  auto pointer_negz = new G4PVPlacement(G4Transform3D(ladderrotation, G4ThreeVector(posx, posy, -m_PosZ[inttlayer][itype])), ladder_volume,
1032  (boost::format("ladder_%d_%d_%d_negz") % inttlayer % itype % icopy).str(), trackerenvelope, false, 0, OverlapCheck());
1033  auto pointer_posz = new G4PVPlacement(G4Transform3D(ladderrotation, G4ThreeVector(posx, posy, +m_PosZ[inttlayer][itype])), ladder_volume,
1034  (boost::format("ladder_%d_%d_%d_posz") % inttlayer % itype % icopy).str(), trackerenvelope, false, 0, OverlapCheck());
1035  if (m_IsActiveMap.find(inttlayer) != m_IsActiveMap.end())
1036  {
1037  m_ActiveVolumeTuple.insert(std::make_pair(pointer_negz, std::make_tuple(inttlayer, itype, icopy, -1)));
1038  m_ActiveVolumeTuple.insert(std::make_pair(pointer_posz, std::make_tuple(inttlayer, itype, icopy, 1)));
1039  }
1040 
1041  // The net effect of the above manipulations for the Z sensitive ladders is that the center of the sensor is at dphi * icopy and at the requested radius
1042  // That us all that the geometry object needs to know, so no changes to that are necessary
1043 
1044  if (itype != 0)
1045  {
1046  // We have added the outer sensor above, now we add the HDI extension tab to the end of the outer sensor HDI
1047  const double posz_ext = (hdi_z_arr[inttlayer][0] + hdi_z) + hdiext_z / 2.;
1048 
1049  new G4PVPlacement(G4Transform3D(ladderrotation, G4ThreeVector(posx, posy, -posz_ext)), ladderext_volume,
1050  (boost::format("ladderext_%d_%d_%d_negz") % inttlayer % itype % icopy).str(), trackerenvelope, false, 0, OverlapCheck());
1051  new G4PVPlacement(G4Transform3D(ladderrotation, G4ThreeVector(posx, posy, +posz_ext)), ladderext_volume,
1052  (boost::format("ladderext_%d_%d_%d_posz") % inttlayer % itype % icopy).str(), trackerenvelope, false, 0, OverlapCheck());
1053  }
1054 
1055  if (Verbosity() > 100)
1056  std::cout << " Ladder copy " << icopy << " radius " << radius << " phi " << phi << " itype " << itype << " posz " << m_PosZ[inttlayer][itype]
1057  << " fRotate " << fRotate << " posx " << posx << " posy " << posy
1058  << std::endl;
1059 
1060  } // end loop over ladder copy placement in phi and positive and negative Z
1061  } // end loop over inner or outer sensor
1062  } // end loop over layers
1063 
1064  // Finally, we add some support material for the silicon detectors
1065 
1066  //
1067  /*
1068  4 rails, which are 12mm OD and 9mm ID tubes at a radius of 168.5 mm. They are spaced equidistantly in phi.
1069  The rails run along the entire length of the TPC and even stick out of the TPC, but I think for the moment you don't have to put the parts that stick out in the simulation.
1070  An inner skin with a ID at 62.416 mm and a thickness of 0.250 mm.
1071  An outer skin with a ID at 120.444 mm and a sandwich of 0.25 mm cfc, 1.5 mm foam and 0.25 mm cfc.
1072 
1073  All of the above are carbon fiber.
1074  */
1076 
1077  // rails
1078  G4Tubs *rail_tube = new G4Tubs((boost::format("si_support_rail")).str(),
1079  supportparams->get_double_param("rail_inner_radius") * cm,
1080  supportparams->get_double_param("rail_outer_radius") * cm,
1081  supportparams->get_double_param("rail_length") * cm / 2.0,
1082  0, 2.0 * M_PI);
1083  G4LogicalVolume *rail_volume = new G4LogicalVolume(rail_tube, GetDetectorMaterial("CFRP_INTT"),
1084  "rail_volume", 0, 0, 0);
1085  if (m_IsSupportActive > 0)
1086  {
1088  }
1089  m_DisplayAction->AddVolume(rail_volume, "Rail");
1090 
1091  double rail_dphi = supportparams->get_double_param("rail_dphi") * deg / rad;
1092  double rail_phi_start = supportparams->get_double_param("rail_phi_start") * deg / rad;
1093  double rail_radius = supportparams->get_double_param("rail_radius") * cm;
1094  for (int i = 0; i < 4; i++)
1095  {
1096  double phi = rail_phi_start + i * rail_dphi;
1097 
1098  // place a copy at each rail phi value
1099  const double posx = rail_radius * cos(phi);
1100  const double posy = rail_radius * sin(phi);
1101 
1102  new G4PVPlacement(0, G4ThreeVector(posx, posy, 0.0), rail_volume,
1103  (boost::format("si_support_rail_%d") % i).str(), trackerenvelope, false, 0, OverlapCheck());
1104  }
1105 
1106  // Outer skin
1107  // G4Tubs *outer_skin_cfcin_tube = new G4Tubs("si_outer_skin_cfcin",
1108  // supportparams->get_double_param("outer_skin_cfcin_inner_radius") * cm,
1109  // supportparams->get_double_param("outer_skin_cfcin_outer_radius") * cm,
1110  // supportparams->get_double_param("outer_skin_cfcin_length") * cm / 2.,
1111  // 0, 2.0 * M_PI);
1112  // G4LogicalVolume *outer_skin_cfcin_volume = new G4LogicalVolume(outer_skin_cfcin_tube, GetDetectorMaterial("CFRP_INTT"),
1113  // "outer_skin_cfcin_volume", 0, 0, 0);
1114 
1115  // G4Tubs *outer_skin_foam_tube = new G4Tubs("si_outer_skin_foam",
1116  // supportparams->get_double_param("outer_skin_foam_inner_radius") * cm,
1117  // supportparams->get_double_param("outer_skin_foam_outer_radius") * cm,
1118  // supportparams->get_double_param("outer_skin_foam_length") * cm / 2.,
1119  // 0, 2.0 * M_PI);
1120  // G4LogicalVolume *outer_skin_foam_volume = new G4LogicalVolume(outer_skin_foam_tube, GetDetectorMaterial("ROHACELL_FOAM_110"),
1121  // "outer_skin_foam_volume", 0, 0, 0);
1122 
1123  // G4Tubs *outer_skin_cfstd::cout_tube = new G4Tubs("si_outer_skin_cfstd::cout",
1124  // supportparams->get_double_param("outer_skin_cfstd::cout_inner_radius") * cm,
1125  // supportparams->get_double_param("outer_skin_cfstd::cout_outer_radius") * cm,
1126  // supportparams->get_double_param("outer_skin_cfstd::cout_length") * cm / 2.,
1127  // 0, 2.0 * M_PI);
1128  // G4LogicalVolume *outer_skin_cfstd::cout_volume = new G4LogicalVolume(outer_skin_cfstd::cout_tube, GetDetectorMaterial("CFRP_INTT"),
1129  // "outer_skin_cfstd::cout_volume", 0, 0, 0);
1130 
1131  G4Tubs *outer_skin_tube = new G4Tubs("si_outer_skin",
1132  supportparams->get_double_param("outer_skin_inner_radius") * cm,
1133  supportparams->get_double_param("outer_skin_outer_radius") * cm,
1134  supportparams->get_double_param("outer_skin_length") * cm / 2.,
1135  0, 2.0 * M_PI);
1136  G4LogicalVolume *outer_skin_volume = new G4LogicalVolume(outer_skin_tube, GetDetectorMaterial("CFRP_INTT"),
1137  "outer_skin_volume", 0, 0, 0);
1138 
1139  if (m_IsSupportActive > 0)
1140  {
1141  // m_PassiveVolumeTuple.insert(std::make_pair(outer_skin_cfcin_volume, std::make_tuple(PHG4InttDefs::SUPPORT_DETID, PHG4InttDefs::INTT_OUTER_SKIN)));
1142  // m_PassiveVolumeTuple.insert(std::make_pair(outer_skin_foam_volume, std::make_tuple(PHG4InttDefs::SUPPORT_DETID, PHG4InttDefs::INTT_OUTER_SKIN)));
1143  // m_PassiveVolumeTuple.insert(std::make_pair(outer_skin_cfstd::cout_volume, std::make_tuple(PHG4InttDefs::SUPPORT_DETID, PHG4InttDefs::INTT_OUTER_SKIN)));
1145  }
1146  // m_DisplayAction->AddVolume(outer_skin_cfcin_volume, "Skin");
1147  // m_DisplayAction->AddVolume(outer_skin_foam_volume, "Skin");
1148  // m_DisplayAction->AddVolume(outer_skin_cfstd::cout_volume, "Skin");
1149  // new G4PVPlacement(0, G4ThreeVector(0, 0.0, 0), outer_skin_cfcin_volume,
1150  // "si_support_outer_skin_cfcin", trackerenvelope, false, 0, OverlapCheck());
1151  // new G4PVPlacement(0, G4ThreeVector(0, 0.0), outer_skin_foam_volume,
1152  // "si_support_outer_skin_foam", trackerenvelope, false, 0, OverlapCheck());
1153  // new G4PVPlacement(0, G4ThreeVector(0, 0.0), outer_skin_cfstd::cout_volume,
1154  // "si_support_outer_skin_cfstd::cout", trackerenvelope, false, 0, OverlapCheck());
1155  m_DisplayAction->AddVolume(outer_skin_volume, "Skin");
1156  new G4PVPlacement(0, G4ThreeVector(0, 0.0, 0), outer_skin_volume,
1157  "si_support_outer_skin_cfcin", trackerenvelope, false, 0, OverlapCheck());
1158 
1159  // Inner skin
1160  G4Tubs *inner_skin_tube = new G4Tubs("si_inner_skin",
1161  supportparams->get_double_param("inner_skin_inner_radius") * cm,
1162  supportparams->get_double_param("inner_skin_outer_radius") * cm,
1163  supportparams->get_double_param("inner_skin_length") * cm / 2.,
1164  0, 2.0 * M_PI);
1165  G4LogicalVolume *inner_skin_volume = new G4LogicalVolume(inner_skin_tube, GetDetectorMaterial("CFRP_INTT"),
1166  "inner_skin_volume", 0, 0, 0);
1167  if (m_IsSupportActive > 0)
1168  {
1170  }
1171  m_DisplayAction->AddVolume(inner_skin_volume, "Skin");
1172 
1173  new G4PVPlacement(0, G4ThreeVector(0, 0.0), inner_skin_volume,
1174  "si_support_inner_skin", trackerenvelope, false, 0, OverlapCheck());
1175 
1176  // Service barrel outer ////////////////////////////////////////////////////////////////////////////////////
1177  G4Tubs *service_barrel_outer_tube = new G4Tubs("si_service_barrel_outer",
1178  supportparams->get_double_param("service_barrel_outer_inner_radius") * cm,
1179  supportparams->get_double_param("service_barrel_outer_outer_radius") * cm,
1180  supportparams->get_double_param("service_barrel_outer_length") * cm / 2.,
1181  0, 2.0 * M_PI);
1182  G4LogicalVolume *service_barrel_outer_volume = new G4LogicalVolume(service_barrel_outer_tube, GetDetectorMaterial("CFRP_INTT"),
1183  "service_barrel_outer_volume", 0, 0, 0);
1184  if (m_IsSupportActive > 0)
1185  {
1186  m_PassiveVolumeTuple.insert(std::make_pair(service_barrel_outer_volume, std::make_tuple(PHG4InttDefs::SUPPORT_DETID, PHG4InttDefs::SERVICE_BARREL_OUTER)));
1187  }
1188  m_DisplayAction->AddVolume(service_barrel_outer_volume, "Skin");
1189 
1190  new G4PVPlacement(0, G4ThreeVector(0, 0.0), service_barrel_outer_volume,
1191  "si_support_service_barrel_outer", trackerenvelope, false, 0, OverlapCheck());
1192 
1193  // Support Tube ////////////////////////////////////////////////////////////////////////////////////
1194  G4Tubs *support_tube_tube = new G4Tubs("si_support_tube",
1195  supportparams->get_double_param("support_tube_inner_radius") * cm,
1196  supportparams->get_double_param("support_tube_outer_radius") * cm,
1197  supportparams->get_double_param("support_tube_length") * cm / 2.,
1198  0, 2.0 * M_PI);
1199  G4LogicalVolume *support_tube_volume = new G4LogicalVolume(support_tube_tube, GetDetectorMaterial("CFRP_INTT"),
1200  "support_tube_volume", 0, 0, 0);
1201  if (m_IsSupportActive > 0)
1202  {
1203  m_PassiveVolumeTuple.insert(std::make_pair(support_tube_volume, std::make_tuple(PHG4InttDefs::SUPPORT_DETID, PHG4InttDefs::SUPPORT_TUBE)));
1204  }
1205  m_DisplayAction->AddVolume(support_tube_volume, "Skin");
1206 
1207  new G4PVPlacement(0, G4ThreeVector(0, 0.0), support_tube_volume,
1208  "si_support_support_tube", trackerenvelope, false, 0, OverlapCheck());
1209 
1210  // Endcap ring in simulations = Endcap rings + endcap staves
1211  int inttlayer = (m_LayerBeginEndIteratorPair.first)->second;
1212  const PHParameters *params1 = m_ParamsContainer->GetParameters(inttlayer);
1213  const int laddertype = params1->get_int_param("laddertype");
1214  const PHParameters *params = m_ParamsContainer->GetParameters(laddertype);
1215 
1216  // Carbon PEEK ring ////////////////////////////////////////////////////////////////////////////////////
1217  G4Tubs *endcap_CP_ring = new G4Tubs("endcap_CP_ring",
1218  supportparams->get_double_param("endcap_CPring_inner_radius") * cm,
1219  supportparams->get_double_param("endcap_CPring_outer_radius") * cm,
1220  supportparams->get_double_param("endcap_CPring_length") * cm / 2.,
1221  0, 2.0 * M_PI);
1222 
1223  G4LogicalVolume *endcap_CP_ring_volume = new G4LogicalVolume(endcap_CP_ring, GetDetectorMaterial("CF30_PEEK70"),
1224  "endcap_CP_ring_volume", 0, 0, 0);
1225  m_DisplayAction->AddVolume(endcap_CP_ring_volume, "EndcapCPRing");
1226 
1227  // new Al-PEEK ring from Jan/2021 //////////////////////////////////////////////////////////////////////////
1228  // outermost part (Al)
1229  G4Tubs *endcap_AlPEEK_Alring_1 = new G4Tubs("endcap_AlPEEK_Alring_1",
1230  supportparams->get_double_param("endcap_AlPEEK_Cring_1_outer_radius") * cm,
1231  supportparams->get_double_param("endcap_AlPEEK_Alring_1_outer_radius") * cm,
1232  supportparams->get_double_param("endcap_AlPEEK_Alring_length") * cm / 2.,
1233  0, 2.0 * M_PI);
1234 
1235  G4LogicalVolume *endcap_AlPEEK_Alring_1_volume = new G4LogicalVolume(endcap_AlPEEK_Alring_1, GetDetectorMaterial("G4_Al"),
1236  "endcap_AlPEEK_Alring_1_volume", 0, 0, 0);
1237  m_DisplayAction->AddVolume(endcap_AlPEEK_Alring_1_volume, "EndcapAlPEEK_Al1");
1238 
1239  // 2nd outermost part (Carbon PEEK)
1240  G4Tubs *endcap_AlPEEK_Cring_1 = new G4Tubs("endcap_AlPEEK_Cring_1",
1241  supportparams->get_double_param("endcap_AlPEEK_Alring_2_outer_radius") * cm,
1242  supportparams->get_double_param("endcap_AlPEEK_Cring_1_outer_radius") * cm,
1243  supportparams->get_double_param("endcap_AlPEEK_Cring_length") * cm / 2.,
1244  0, 2.0 * M_PI);
1245 
1246  G4LogicalVolume *endcap_AlPEEK_Cring_1_volume = new G4LogicalVolume(endcap_AlPEEK_Cring_1, GetDetectorMaterial("CF30_PEEK70"),
1247  "endcap_AlPEEK_Cring_1_volume", 0, 0, 0);
1248  m_DisplayAction->AddVolume(endcap_AlPEEK_Cring_1_volume, "EndcapAlPEEK_C1");
1249 
1250  // 3rd outermost part (Al)
1251  G4Tubs *endcap_AlPEEK_Alring_2 = new G4Tubs("endcap_AlPEEK_Alring_2",
1252  supportparams->get_double_param("endcap_AlPEEK_Cring_2_outer_radius") * cm,
1253  supportparams->get_double_param("endcap_AlPEEK_Alring_2_outer_radius") * cm,
1254  supportparams->get_double_param("endcap_AlPEEK_Alring_length") * cm / 2.,
1255  0, 2.0 * M_PI);
1256 
1257  G4LogicalVolume *endcap_AlPEEK_Alring_2_volume = new G4LogicalVolume(endcap_AlPEEK_Alring_2, GetDetectorMaterial("G4_Al"),
1258  "endcap_AlPEEK_Alring_2_volume", 0, 0, 0);
1259  m_DisplayAction->AddVolume(endcap_AlPEEK_Alring_2_volume, "EndcapAlPEEK_Al2");
1260 
1261  // 4th outermost part (Carbon PEEK)
1262  G4Tubs *endcap_AlPEEK_Cring_2 = new G4Tubs("endcap_AlPEEK_Cring_2",
1263  supportparams->get_double_param("endcap_AlPEEK_Alring_3_outer_radius") * cm,
1264  supportparams->get_double_param("endcap_AlPEEK_Cring_2_outer_radius") * cm,
1265  supportparams->get_double_param("endcap_AlPEEK_Cring_length") * cm / 2.,
1266  0, 2.0 * M_PI);
1267 
1268  G4LogicalVolume *endcap_AlPEEK_Cring_2_volume = new G4LogicalVolume(endcap_AlPEEK_Cring_2, GetDetectorMaterial("CF30_PEEK70"),
1269  "endcap_AlPEEK_Cring_2_volume", 0, 0, 0);
1270  m_DisplayAction->AddVolume(endcap_AlPEEK_Cring_2_volume, "EndcapAlPEEK_C2");
1271 
1272  // 5th outermost part (innermost) (Al)
1273  G4Tubs *endcap_AlPEEK_Alring_3 = new G4Tubs("endcap_AlPEEK_Alring_3",
1274  supportparams->get_double_param("endcap_AlPEEK_Alring_3_inner_radius") * cm,
1275  supportparams->get_double_param("endcap_AlPEEK_Alring_3_outer_radius") * cm,
1276  supportparams->get_double_param("endcap_AlPEEK_Alring_length") * cm / 2.,
1277  0, 2.0 * M_PI);
1278 
1279  G4LogicalVolume *endcap_AlPEEK_Alring_3_volume = new G4LogicalVolume(endcap_AlPEEK_Alring_3, GetDetectorMaterial("G4_Al"),
1280  "endcap_AlPEEK_Alring_3_volume", 0, 0, 0);
1281  m_DisplayAction->AddVolume(endcap_AlPEEK_Alring_3_volume, "EndcapAlPEEK_Al3");
1282 
1283  if (m_IsEndcapActive)
1284  {
1285  double endcap_outer_edge_z = 0.0; // absolute z-coordinate of outer edge (furthest side from the origin) of the endcap, used for bus extender
1286  if (supportparams->get_int_param("endcap_ring_type") == 0) // Place Al endcap rings
1287  {
1288  // Aluminum ring
1289  G4Tubs *endcap_Al_ring = new G4Tubs("endcap_Al_ring",
1290  supportparams->get_double_param("endcap_Alring_inner_radius") * cm,
1291  supportparams->get_double_param("endcap_Alring_outer_radius") * cm,
1292  supportparams->get_double_param("endcap_Alring_length") * cm / 2.,
1293  0, 2.0 * M_PI);
1294 
1295  G4LogicalVolume *endcap_Al_ring_volume = new G4LogicalVolume(endcap_Al_ring, GetDetectorMaterial("Al6061T6"),
1296  "endcap_Al_ring_volume", 0, 0, 0);
1297 
1298  // Stainlees steal ring
1299  G4Tubs *endcap_SS_ring = new G4Tubs("endcap_SS_ring",
1300  supportparams->get_double_param("endcap_SSring_inner_radius") * cm,
1301  supportparams->get_double_param("endcap_SSring_outer_radius") * cm,
1302  supportparams->get_double_param("endcap_SSring_length") * cm / 2.,
1303  0, 2.0 * M_PI);
1304 
1305  G4LogicalVolume *endcap_SS_ring_volume = new G4LogicalVolume(endcap_SS_ring, GetDetectorMaterial("SS316"),
1306  "endcap_SS_ring_volume", 0, 0, 0);
1307 
1308  // Water Glycol ring
1309  G4Tubs *endcap_WG_ring = new G4Tubs("endcap_WG_ring",
1310  supportparams->get_double_param("endcap_WGring_inner_radius") * cm,
1311  supportparams->get_double_param("endcap_WGring_outer_radius") * cm,
1312  supportparams->get_double_param("endcap_WGring_length") * cm / 2.,
1313  0, 2.0 * M_PI);
1314 
1315  G4LogicalVolume *endcap_WG_ring_volume = new G4LogicalVolume(endcap_WG_ring, GetDetectorMaterial("WaterGlycol_INTT"),
1316  "endcap_WG_ring_volume", 0, 0, 0);
1317 
1318  double endcap_ring_z = supportparams->get_double_param("endcap_ring_z") * cm;
1319  for (int i = 0; i < 2; i++) // i=0 : positive z, i=1 negative z
1320  {
1321  endcap_ring_z = (i == 0) ? endcap_ring_z : -1.0 * endcap_ring_z;
1322 
1323  double width_WGring_z = supportparams->get_double_param("endcap_WGring_length") * cm;
1324  double width_SSring_z = supportparams->get_double_param("endcap_SSring_length") * cm;
1325  double width_Alring_z = supportparams->get_double_param("endcap_Alring_length") * cm;
1326 
1327  for (int j = 0; j < 2; j++) // j=0 : positive side z, j=1 negative side z
1328  {
1329  width_WGring_z = (j == 0) ? width_WGring_z : -1.0 * width_WGring_z;
1330  width_SSring_z = (j == 0) ? width_SSring_z : -1.0 * width_SSring_z;
1331  width_Alring_z = (j == 0) ? width_Alring_z : -1.0 * width_Alring_z;
1332 
1333  double cent_WGring_z = endcap_ring_z + width_WGring_z / 2.;
1334  double cent_SSring_z = endcap_ring_z + width_WGring_z + width_SSring_z / 2.;
1335  double cent_Alring_z = endcap_ring_z + width_WGring_z + width_SSring_z + width_Alring_z / 2.;
1336  endcap_outer_edge_z = fabs(endcap_ring_z) + fabs(width_WGring_z + width_SSring_z + width_Alring_z / 2.); // absolute distance from origin
1337 
1338  new G4PVPlacement(0, G4ThreeVector(0, 0, cent_WGring_z),
1339  endcap_WG_ring_volume,
1340  (boost::format("endcap_WG_ring_pv_%d_%d") % i % j).str(),
1341  trackerenvelope, false, 0, OverlapCheck());
1342 
1343  new G4PVPlacement(0, G4ThreeVector(0, 0, cent_SSring_z),
1344  endcap_SS_ring_volume,
1345  (boost::format("endcap_SS_ring_pv_%d_%d") % i % j).str(),
1346  trackerenvelope, false, 0, OverlapCheck());
1347 
1348  new G4PVPlacement(0, G4ThreeVector(0, 0, cent_Alring_z),
1349  endcap_Al_ring_volume,
1350  (boost::format("endcap_Al_ring_pv_%d_%d") % i % j).str(),
1351  trackerenvelope, false, 0, OverlapCheck());
1352  } // end of the loop over positive/negative sides with j
1353  } // end of the loop over positive/negative sides with i
1354  } // end of endcap_ring_type == 0
1355  else if (supportparams->get_int_param("endcap_ring_type") == 1) // Place CP endcap rings
1356  {
1357  double endcap_ring_z = supportparams->get_double_param("endcap_CPring_z") * cm;
1358 
1359  for (int i = 0; i < 2; i++) // i=0 : positive z, i=1 negative z
1360  {
1361  endcap_ring_z = (i == 0) ? endcap_ring_z : -1.0 * endcap_ring_z;
1362  endcap_outer_edge_z = fabs(endcap_ring_z);
1363 
1364  new G4PVPlacement(0, G4ThreeVector(0, 0, endcap_ring_z),
1365  endcap_CP_ring_volume,
1366  (boost::format("endcap_CP_ring_pv_%d") % i).str(),
1367  trackerenvelope, false, 0, OverlapCheck());
1368 
1369  } // end of the loop over positive/negative sides
1370  } // end of endcap_ring_type == 1
1371  else if (supportparams->get_int_param("endcap_ring_type") == 2) // the new endcap introduced in Jan/2021
1372  {
1373  double si_0_width = params->get_double_param("strip_z_0") * params->get_int_param("nstrips_z_sensor_0") * cm + 2 * params->get_double_param("sensor_edge_z") * cm; // length of the smaller cells
1374  double si_1_width = params->get_double_param("strip_z_1") * params->get_int_param("nstrips_z_sensor_1") * cm + 2 * params->get_double_param("sensor_edge_z") * cm; // length of the larger cells
1375  double sifull_width = si_0_width + si_1_width; // length of the Si module
1376  double hdi_width = sifull_width + params->get_double_param("hdi_edge_z") * cm;
1377  double hdiext_width = params->get_double_param("halfladder_inside_z") * cm - sifull_width;
1378  double hdifull_width = hdi_width + hdiext_width; // length of a half lader
1379  double endcap_ring_z = hdifull_width + supportparams->get_double_param("endcap_AlPEEK_Cring_length") / 2.0 * cm;
1380 
1381  for (int i = 0; i < 2; i++) // i=0 : positive z, i=1 negative z
1382  {
1383  endcap_ring_z = (i == 0) ? endcap_ring_z : -1.0 * endcap_ring_z;
1384  endcap_outer_edge_z = fabs(endcap_ring_z) + supportparams->get_double_param("endcap_AlPEEK_Alring_length") * cm / 2.0;
1385 
1386  new G4PVPlacement(0, G4ThreeVector(0, 0, endcap_ring_z),
1387  endcap_AlPEEK_Alring_1_volume,
1388  (boost::format("endcap_AlPEEK_Alring_1_pv_%d") % i).str(),
1389  trackerenvelope, false, 0, OverlapCheck());
1390 
1391  new G4PVPlacement(0, G4ThreeVector(0, 0, endcap_ring_z),
1392  endcap_AlPEEK_Cring_1_volume,
1393  (boost::format("endcap_AlPEEK_Cring_1_pv_%d") % i).str(),
1394  trackerenvelope, false, 0, OverlapCheck());
1395 
1396  new G4PVPlacement(0, G4ThreeVector(0, 0, endcap_ring_z),
1397  endcap_AlPEEK_Alring_2_volume,
1398  (boost::format("endcap_AlPEEK_Alring_2_pv_%d") % i).str(),
1399  trackerenvelope, false, 0, OverlapCheck());
1400 
1401  new G4PVPlacement(0, G4ThreeVector(0, 0, endcap_ring_z),
1402  endcap_AlPEEK_Cring_2_volume,
1403  (boost::format("endcap_AlPEEK_Cring_2_pv_%d") % i).str(),
1404  trackerenvelope, false, 0, OverlapCheck());
1405 
1406  new G4PVPlacement(0, G4ThreeVector(0, 0, endcap_ring_z),
1407  endcap_AlPEEK_Alring_3_volume,
1408  (boost::format("endcap_AlPEEK_Alring_3_pv_%d") % i).str(),
1409  trackerenvelope, false, 0, OverlapCheck());
1410  }
1411  }
1412 
1414  // Mimic cylinders for the bus extender (very simplified for the moment)
1415  double bus_extender_radius_innermost = supportparams->get_double_param("bus_extender_radius") * cm;
1416  double bus_extender_length = supportparams->get_double_param("bus_extender_length") * cm;
1417  double bus_extender_copper_x = supportparams->get_double_param("bus_extender_copper_x") * cm;
1418  double bus_extender_kapton_x = supportparams->get_double_param("bus_extender_kapton_x") * cm;
1419 
1420  // copper layer of the bus extender for the inner barrel
1421  double inner_radius = bus_extender_radius_innermost;
1422  G4Tubs *bus_extender_copper_inner = new G4Tubs("bus_extender_coppe_innerr",
1423  inner_radius, inner_radius + bus_extender_copper_x,
1424  bus_extender_length / 2.0, 0, 2.0 * M_PI);
1425 
1426  G4LogicalVolume *bus_extender_copper_inner_volume = new G4LogicalVolume(bus_extender_copper_inner, GetDetectorMaterial("G4_Cu"),
1427  "bus_extender_copper_inner_volume", 0, 0, 0);
1428  m_DisplayAction->AddVolume(bus_extender_copper_inner_volume, "BusExtenderCopperInner");
1429 
1430  // kapton layer of the bus extender for the inner barrel
1431  inner_radius += bus_extender_copper_x;
1432  G4Tubs *bus_extender_kapton_inner = new G4Tubs("bus_extender_kapton_inner",
1433  inner_radius, inner_radius + bus_extender_kapton_x,
1434  bus_extender_length / 2.0, 0, 2.0 * M_PI);
1435 
1436  G4LogicalVolume *bus_extender_kapton_inner_volume = new G4LogicalVolume(bus_extender_kapton_inner, GetDetectorMaterial("G4_KAPTON"),
1437  "bus_extender_kapton_inner_volume", 0, 0, 0);
1438  m_DisplayAction->AddVolume(bus_extender_kapton_inner_volume, "BusExtenderKaptonInner");
1439 
1440  // copper layer of the bus extender for the outer outerbarrel
1441  inner_radius += bus_extender_kapton_x;
1442  G4Tubs *bus_extender_copper_outer = new G4Tubs("bus_extender_copper_outer",
1443  inner_radius, inner_radius + bus_extender_copper_x,
1444  bus_extender_length / 2.0, 0, 2.0 * M_PI);
1445 
1446  G4LogicalVolume *bus_extender_copper_outer_volume = new G4LogicalVolume(bus_extender_copper_outer, GetDetectorMaterial("G4_Cu"),
1447  "bus_extender_copper_outer_volume", 0, 0, 0);
1448  m_DisplayAction->AddVolume(bus_extender_copper_outer_volume, "BusExtenderCopperOuter");
1449 
1450  // kapton layer of the bus extender for the outer barrel
1451  inner_radius += bus_extender_copper_x;
1452  G4Tubs *bus_extender_kapton_outer = new G4Tubs("bus_extender_kapton_outer",
1453  inner_radius, inner_radius + bus_extender_kapton_x,
1454  bus_extender_length / 2.0, 0, 2.0 * M_PI);
1455 
1456  G4LogicalVolume *bus_extender_kapton_outer_volume = new G4LogicalVolume(bus_extender_kapton_outer, GetDetectorMaterial("G4_KAPTON"),
1457  "bus_extender_kapton_outer_volume", 0, 0, 0);
1458  m_DisplayAction->AddVolume(bus_extender_kapton_outer_volume, "BusExtenderKaptonOuter");
1459 
1460  double bus_extender_z = endcap_outer_edge_z;
1461  for (int i = 0; i < 2; i++) // loop over both positive and negative sides to put the cylinders, i=0 : positive z, i=1 negative z
1462  {
1463  // copper layer of bus extender for the inner barrel
1464  double cent_bus_extender_z = bus_extender_z + bus_extender_length / 2.0;
1465  cent_bus_extender_z *= (i == 0) ? 1.0 : -1.0;
1466  new G4PVPlacement(0, G4ThreeVector(0, 0, cent_bus_extender_z),
1467  bus_extender_copper_inner_volume,
1468  (boost::format("bus_extender_copper_inner_layer_pv_%d") % i).str(),
1469  trackerenvelope, false, 0, OverlapCheck());
1470 
1471  // kapton layer of bus extender for the inner barrel
1472  new G4PVPlacement(0, G4ThreeVector(0, 0, cent_bus_extender_z),
1473  bus_extender_kapton_inner_volume,
1474  (boost::format("bus_extender_kapton_inner_layer_pv_%d") % i).str(),
1475  trackerenvelope, false, 0, OverlapCheck());
1476 
1477  // copper layer of bus extender for the outer barrel
1478  new G4PVPlacement(0, G4ThreeVector(0, 0, cent_bus_extender_z),
1479  bus_extender_copper_outer_volume,
1480  (boost::format("bus_extender_copper_outer_layer_pv_%d") % i).str(),
1481  trackerenvelope, false, 0, OverlapCheck());
1482 
1483  // kapton layer of bus extender for the outer barrel
1484  new G4PVPlacement(0, G4ThreeVector(0, 0, cent_bus_extender_z),
1485  bus_extender_kapton_outer_volume,
1486  (boost::format("bus_extender_kapton_outer_layer_pv_%d") % i).str(),
1487  trackerenvelope, false, 0, OverlapCheck());
1488  }
1489  }
1490 
1491  return 0;
1492 }
1493 
1495 {
1496  int active = 0;
1497  // std::map<int, int>::const_iterator iter;
1498  for (auto iter = m_IsActiveMap.begin(); iter != m_IsActiveMap.end(); ++iter)
1499  {
1500  if (iter->second > 0)
1501  {
1502  active = 1;
1503  break;
1504  }
1505  }
1506  if (active)
1507  {
1508  std::string geonode = (m_SuperDetector != "NONE") ? (boost::format("CYLINDERGEOM_%s") % m_SuperDetector).str() : (boost::format("CYLINDERGEOM_%s") % m_DetectorType).str();
1509 
1510  PHG4CylinderGeomContainer *geo = findNode::getClass<PHG4CylinderGeomContainer>(topNode(), geonode);
1511  if (!geo)
1512  {
1513  geo = new PHG4CylinderGeomContainer();
1514  PHNodeIterator iter(topNode());
1515  PHCompositeNode *runNode = dynamic_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "RUN"));
1516  PHIODataNode<PHObject> *newNode = new PHIODataNode<PHObject>(geo, geonode, "PHObject");
1517  runNode->addNode(newNode);
1518  }
1519 
1520  for (auto layeriter = m_LayerBeginEndIteratorPair.first; layeriter != m_LayerBeginEndIteratorPair.second; ++layeriter)
1521  {
1522  const int sphxlayer = layeriter->first;
1523  const int inttlayer = layeriter->second;
1524  int ilayer = inttlayer;
1525  const PHParameters *params_layer = m_ParamsContainer->GetParameters(inttlayer);
1526  const int laddertype = params_layer->get_int_param("laddertype");
1527  // parameters are stored in cm per our convention
1528  const PHParameters *params = m_ParamsContainer->GetParameters(laddertype);
1529  CylinderGeomIntt *mygeom = new CylinderGeomIntt(
1530  sphxlayer,
1531  params->get_double_param("strip_x"),
1532  params->get_double_param("strip_y"),
1533  params->get_double_param("strip_z_0"),
1534  params->get_double_param("strip_z_1"),
1535  params->get_int_param("nstrips_z_sensor_0"),
1536  params->get_int_param("nstrips_z_sensor_1"),
1537  params->get_int_param("nstrips_phi_sensor"),
1538  params_layer->get_int_param("nladder"),
1539  m_PosZ[ilayer][0] / cm, // m_PosZ uses G4 internal units, needs to be converted to cm
1540  m_PosZ[ilayer][1] / cm,
1541  m_SensorRadius[ilayer] / cm,
1542  0.0,
1543  params_layer->get_double_param("offsetphi") * deg / rad, // expects radians
1544  params_layer->get_double_param("offsetrot") * deg / rad);
1545  geo->AddLayerGeom(sphxlayer, mygeom);
1546  if (Verbosity() > 0)
1547  {
1548  geo->identify();
1549  }
1550  }
1551  }
1552 }
1553 
1554 std::map<G4VPhysicalVolume *, std::tuple<int, int, int, int>>::const_iterator
1556 {
1557  auto iter = m_ActiveVolumeTuple.find(physvol);
1558  if (iter == m_ActiveVolumeTuple.end())
1559  {
1560  std::cout << PHWHERE << " Volume " << physvol->GetName() << " not in active volume tuple" << std::endl;
1561  gSystem->Exit(1);
1562  }
1563  return iter;
1564 }
1565 
1566 std::map<G4LogicalVolume *, std::tuple<int, int>>::const_iterator
1567 PHG4InttDetector::get_PassiveVolumeTuple(G4LogicalVolume *logvol) const
1568 {
1569  auto iter = m_PassiveVolumeTuple.find(logvol);
1570  if (iter == m_PassiveVolumeTuple.end())
1571  {
1572  std::cout << PHWHERE << " Volume " << logvol->GetName() << " not in passive volume tuple" << std::endl;
1573  gSystem->Exit(1);
1574  }
1575  return iter;
1576 }