ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
PHG4FCalDetector.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file PHG4FCalDetector.cc
1 #include "PHG4FCalDetector.h"
2 
4 
5 #include <g4main/PHG4Detector.h> // for PHG4Detector
6 
7 #include <Geant4/G4Box.hh>
8 #include <Geant4/G4Colour.hh>
9 #include <Geant4/G4LogicalVolume.hh>
10 #include <Geant4/G4Material.hh>
11 #include <Geant4/G4NistManager.hh>
12 #include <Geant4/G4PVPlacement.hh>
13 #include <Geant4/G4Region.hh> // for G4Region
14 #include <Geant4/G4String.hh> // for G4String
15 #include <Geant4/G4SystemOfUnits.hh>
16 #include <Geant4/G4ThreeVector.hh> // for G4ThreeVector
17 #include <Geant4/G4Types.hh>
18 #include <Geant4/G4VisAttributes.hh>
19 
20 #include <cstdlib> // for NULL, exit
21 #include <iostream> // for stringstream, operator<<
22 #include <map>
23 #include <sstream>
24 #include <utility> // for pair
25 
26 class PHCompositeNode;
27 
28 using namespace std;
29 
31  : PHG4Detector(subsys, Node, name)
32  , length(1.0 * m)
33  , absorber_thickness(5.0 * mm)
34  , scintillator_thickness(5.0 * cm)
35  , nlayers(4)
36  , segments_per_column(20)
37  , segments_per_thickness(1)
38  , z_position(100.0 * cm)
39  , layer_separation(1.0 * mm)
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 PHG4FCalDetector::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 / ((double) segments_per_thickness);
61  unsigned int z_index = (unsigned int) ((z - z_start) / z_width);
62 
63  double xy_start = -0.5 * length;
64  double xy_width = length / ((double) segments_per_column);
65 
66  unsigned int x_index = (unsigned int) ((x - xy_start) / xy_width);
67  unsigned int y_index = (unsigned int) ((y - xy_start) / xy_width);
68 
69  zcenter = z_start + z_index * z_width + 0.5 * z_width;
70  xcenter = xy_start + x_index * xy_width + 0.5 * xy_width;
71  ycenter = xy_start + y_index * xy_width + 0.5 * xy_width;
72 
73  return (z_index * (segments_per_column * segments_per_column) + y_index * segments_per_column + x_index);
74 }
75 
77 {
78  // const G4MaterialTable* mattab = GetDetectorMaterialTable();
79  // for(unsigned int i=0;i<mattab->size();i++)
80  // {
81  // cout<<mattab->at(i)->GetName()<<endl;
82  // }
83 
85  ScintillatorMaterial = man->FindOrBuildMaterial("G4_POLYSTYRENE");
86  AbsorberMaterial = man->FindOrBuildMaterial("G4_Pb");
87 
88  if (!ScintillatorMaterial ||
90  {
91  std::cout << "PHG4SvxDetector::Construct - Error: Can not set material" << std::endl;
92  exit(-1);
93  }
94 
95  G4VisAttributes* polystyreneVis = new G4VisAttributes(G4Colour(0.0, 1.0, 1.0));
96  polystyreneVis->SetVisibility(true);
97  polystyreneVis->SetForceSolid(true);
98 
99  G4VisAttributes* leadVis = new G4VisAttributes(G4Colour(0.925, 0.345, 0));
100  leadVis->SetVisibility(true);
101  leadVis->SetForceSolid(true);
102 
103  _region = new G4Region("FCALREGION");
105 
106  G4double current_center_position = z_position;
107  for (unsigned int layer = 0; layer < nlayers; layer++)
108  {
109  stringstream ss;
110 
111  ss.clear();
112  ss.str("");
113  ss << "FCalAbsSolid_" << layer;
114  absorber_solid_[layer] = new G4Box(G4String(ss.str()), 0.5 * length, 0.5 * length, 0.5 * absorber_thickness);
115 
116  ss.clear();
117  ss.str("");
118  ss << "FCalAbsLogical_" << layer;
120  absorber_logic_[layer]->SetVisAttributes(leadVis);
122 
123  ss.clear();
124  ss.str("");
125  ss << "FCalAbsPhysical_" << layer;
126 
127  absorber_physi_[layer] = new G4PVPlacement(0, G4ThreeVector(0. * cm, 0. * cm, current_center_position), absorber_logic_[layer], G4String(ss.str()), logicWorld, false, 0, false);
128 
129  current_center_position += 0.5 * absorber_thickness;
130  current_center_position += 0.5 * scintillator_thickness;
131  current_center_position += layer_separation;
132 
133  ss.clear();
134  ss.str("");
135  ss << "FCalScintSolid_" << layer;
136  scintillator_solid_[layer] = new G4Box(G4String(ss.str()), 0.5 * length, 0.5 * length, 0.5 * scintillator_thickness);
137 
138  ss.clear();
139  ss.str("");
140  ss << "FCalScintLogical_" << layer;
142  scintillator_logic_[layer]->SetVisAttributes(polystyreneVis);
144 
145  ss.clear();
146  ss.str("");
147  ss << "FCalScintPhysical_" << layer;
148 
149  scintillator_physi_[layer] = new G4PVPlacement(0, G4ThreeVector(0. * cm, 0. * cm, current_center_position), scintillator_logic_[layer], G4String(ss.str()), logicWorld, false, 0, false);
150 
151  current_center_position += 0.5 * absorber_thickness;
152  current_center_position += 0.5 * scintillator_thickness;
153  current_center_position += layer_separation;
154  }
155 }
156 
158 {
159  //loop over the physical volumes and see if this is a match
160  std::map<unsigned int, G4VPhysicalVolume*>::iterator vol_iter = scintillator_physi_.begin();
161  for (; vol_iter != scintillator_physi_.end(); ++vol_iter)
162  {
163  if (vol_iter->second == volume)
164  return true;
165  }
166  return false;
167 }
168 
170 {
171  //loop over the physical volumes and see if this is a match
172  std::map<unsigned int, G4VPhysicalVolume*>::iterator vol_iter = scintillator_physi_.begin();
173  for (; vol_iter != scintillator_physi_.end(); ++vol_iter)
174  {
175  if (vol_iter->second == volume)
176  return vol_iter->first;
177  }
178  return -1;
179 }