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