ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PHG4FPbScDetector.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file PHG4FPbScDetector.cc
1 #include "PHG4FPbScDetector.h"
2 
3 #include <Geant4/G4Box.hh>
4 #include <Geant4/G4Colour.hh>
5 #include <Geant4/G4LogicalVolume.hh>
6 #include <Geant4/G4Material.hh>
7 #include <Geant4/G4NistManager.hh>
8 #include <Geant4/G4PVPlacement.hh>
9 #include <Geant4/G4Region.hh> // for G4Region
10 #include <Geant4/G4String.hh> // for G4String
11 #include <Geant4/G4SystemOfUnits.hh> // for cm
12 #include <Geant4/G4ThreeVector.hh> // for G4ThreeVector
13 #include <Geant4/G4Types.hh>
14 #include <Geant4/G4VisAttributes.hh>
15 
16 #include <cstdlib> // for NULL, exit
17 #include <iostream> // for stringstream, operator<<
18 #include <map>
19 #include <sstream>
20 #include <utility> // for pair
21 
22 class PHCompositeNode;
23 
24 using namespace std;
25 
27  : PHG4Detector(subsys, Node, nam)
28  , tower_cross_section(5.535 * cm)
29  , segments_per_column(12 * 6)
30  , segments_per_height(12 * 3)
31  , length(tower_cross_section * segments_per_column)
32  , height(tower_cross_section * segments_per_height)
33  , absorber_thickness(0.15 * cm)
34  , scintillator_thickness(0.4 * cm)
35  , nlayers(66)
36  , x_position(0.0 * cm)
37  , y_position(0.0 * cm)
38  , z_position(400.0 * cm)
39  , layer_separation(0.0)
40  , AbsorberMaterial(nullptr)
41  , ScintillatorMaterial(nullptr)
42  , _region(nullptr)
43 {
44 }
45 
47 {
48  // search the material by its name and assign if found
49  if (G4Material* material_ptr = GetDetectorMaterial(material))
50  {
51  return material_ptr;
52  }
53  else
54  return 0;
55 }
56 
57 unsigned int PHG4FPbScDetector::computeIndex(unsigned int layer, G4double x, G4double y, G4double z, G4double& xcenter, G4double& ycenter, G4double& zcenter)
58 {
59  double z_start = z_position + absorber_thickness * (((double) layer) + 1.) + layer_separation * (((double) layer) + 1.) + scintillator_thickness * (((double) layer));
60  double z_width = scintillator_thickness;
61  unsigned int z_index = (unsigned int) ((z - z_start) / z_width);
62 
63  double x_start = x_position - 0.5 * length;
64  double y_start = y_position - 0.5 * height;
65  double xy_width = tower_cross_section;
66 
67  unsigned int x_index = (unsigned int) ((x - x_start) / xy_width);
68  unsigned int y_index = (unsigned int) ((y - y_start) / xy_width);
69 
70  zcenter = z; //z_start + z_index*z_width + 0.5*z_width;
71  xcenter = x; //x_start + x_index*xy_width + 0.5*xy_width;
72  ycenter = y; //y_start + y_index*xy_width + 0.5*xy_width;
73 
74  return (z_index * (segments_per_column * segments_per_height) + y_index * segments_per_column + x_index);
75 }
76 
78 {
79  // const G4MaterialTable* mattab = G4Material::GetMaterialTable();
80  // for(unsigned int i=0;i<mattab->size();i++)
81  // {
82  // cout<<mattab->at(i)->GetName()<<endl;
83  // }
84 
86  ScintillatorMaterial = man->FindOrBuildMaterial("G4_POLYSTYRENE");
87  AbsorberMaterial = man->FindOrBuildMaterial("G4_Pb");
88 
89  if (!ScintillatorMaterial ||
91  {
92  std::cout << "PHG4FPbScDetector::Construct - Error: Can not set material" << std::endl;
93  exit(-1);
94  }
95 
96  G4VisAttributes* polystyreneVis = new G4VisAttributes(G4Colour::Blue());
97  polystyreneVis->SetVisibility(true);
98  polystyreneVis->SetForceSolid(true);
99 
100  G4VisAttributes* leadVis = new G4VisAttributes(G4Colour(0.925, 0.345, 0));
101  leadVis->SetVisibility(true);
102  leadVis->SetForceSolid(true);
103 
104  _region = new G4Region("FPBSCREGION");
105  // _region->SetRegionalSteppingAction(new PHG4FPbScSteppingAction(this));
106 
107  G4double current_center_position = z_position;
108  for (unsigned int layer = 0; layer < nlayers; layer++)
109  {
110  stringstream ss;
111 
112  ss.clear();
113  ss.str("");
114  ss << "FPbScAbsSolid_" << layer;
115  absorber_solid_[layer] = new G4Box(G4String(ss.str()),
116  0.5 * length,
117  0.5 * height,
118  0.5 * absorber_thickness);
119 
120  ss.clear();
121  ss.str("");
122  ss << "FPbScAbsLogical_" << layer;
125  G4String(ss.str()), 0, 0, 0);
126  absorber_logic_[layer]->SetVisAttributes(leadVis);
127  // _region->AddRootLogicalVolume(absorber_logic_[layer]);
128 
129  ss.clear();
130  ss.str("");
131  ss << "FPbScAbsPhysical_" << layer;
132 
135  current_center_position),
136  absorber_logic_[layer],
137  G4String(ss.str()),
138  logicWorld, false, 0, false);
139 
140  current_center_position += 0.5 * absorber_thickness;
141  current_center_position += 0.5 * scintillator_thickness;
142 
143  ss.clear();
144  ss.str("");
145  ss << "FPbScScintSolid_" << layer;
146  scintillator_solid_[layer] = new G4Box(G4String(ss.str()),
147  0.5 * length, 0.5 * height,
148  0.5 * scintillator_thickness);
149 
150  ss.clear();
151  ss.str("");
152  ss << "FPbScScintLogical_" << layer;
154  ScintillatorMaterial, G4String(ss.str()), 0, 0, 0);
155  scintillator_logic_[layer]->SetVisAttributes(polystyreneVis);
156  // _region->AddRootLogicalVolume(scintillator_logic_[layer]);
157 
158  ss.clear();
159  ss.str("");
160  ss << "FPbScScintPhysical_" << layer;
161 
162  scintillator_physi_[layer] = new G4PVPlacement(0, G4ThreeVector(x_position, y_position, current_center_position),
163  scintillator_logic_[layer],
164  G4String(ss.str()),
165  logicWorld, false, 0, false);
166 
167  current_center_position += 0.5 * scintillator_thickness;
168  current_center_position += layer_separation;
169  current_center_position += 0.5 * absorber_thickness;
170  }
171 }
172 
174 {
175  //loop over the physical volumes and see if this is a match
176  std::map<unsigned int, G4VPhysicalVolume*>::iterator vol_iter = scintillator_physi_.begin();
177  for (; vol_iter != scintillator_physi_.end(); ++vol_iter)
178  {
179  if (vol_iter->second == volume)
180  {
181  return true;
182  }
183  }
184  return false;
185 }
186 
188 {
189  //loop over the physical volumes and see if this is a match
190  std::map<unsigned int, G4VPhysicalVolume*>::iterator vol_iter = scintillator_physi_.begin();
191  for (; vol_iter != scintillator_physi_.end(); ++vol_iter)
192  {
193  if (vol_iter->second == volume)
194  return vol_iter->first;
195  }
196  return -1;
197 }