ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AllSiliconTrackerDetector.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file AllSiliconTrackerDetector.cc
2 
5 
6 #include <phparameter/PHParameters.h>
7 
8 #include <g4main/PHG4Detector.h>
10 #include <g4main/PHG4Subsystem.h>
11 
12 #include <phool/PHCompositeNode.h>
13 #include <phool/PHIODataNode.h>
14 #include <phool/PHNode.h>
15 #include <phool/PHNodeIterator.h>
16 #include <phool/PHObject.h>
17 #include <phool/getClass.h>
18 
19 #include <TSystem.h>
20 
21 #include <Geant4/G4AssemblyVolume.hh>
22 #include <Geant4/G4Color.hh>
23 #include <Geant4/G4GDMLParser.hh>
24 #include <Geant4/G4GDMLReadStructure.hh> // for G4GDMLReadStructure
25 #include <Geant4/G4LogicalVolume.hh>
26 #include <Geant4/G4LogicalVolumeStore.hh>
27 #include <Geant4/G4Material.hh>
28 #include <Geant4/G4PVPlacement.hh>
29 #include <Geant4/G4SolidStore.hh>
30 #include <Geant4/G4SystemOfUnits.hh>
31 #include <Geant4/G4VisAttributes.hh>
32 
33 #include <boost/algorithm/string.hpp>
34 #include <boost/lexical_cast.hpp>
35 
36 #include <cmath>
37 #include <iostream>
38 #include <memory>
39 
40 class G4VSolid;
41 class PHCompositeNode;
42 
43 using namespace std;
44 
45 //____________________________________________________________________________..
47  PHCompositeNode *Node,
49  const std::string &dnam)
50  : PHG4Detector(subsys, Node, dnam)
51  , m_DisplayAction(dynamic_cast<AllSiliconTrackerDisplayAction *>(subsys->GetDisplayAction()))
52  , m_Params(parameters)
53  , m_GDMPath(parameters->get_string_param("GDMPath"))
54  , m_Active(m_Params->get_int_param("active"))
55  , m_AbsorberActive(parameters->get_int_param("absorberactive"))
56 {
57 }
58 
59 //_______________________________________________________________
61 {
62  if (m_Active)
63  {
64  auto iter = m_ActivePhysicalVolumesSet.find(volume);
65  if (iter != m_ActivePhysicalVolumesSet.end())
66  {
67  return iter->second;
68  }
69  }
70  if (m_AbsorberActive)
71  {
72  auto iter = m_PassivePhysicalVolumesSet.find(volume);
73  if (iter != m_PassivePhysicalVolumesSet.end())
74  {
75  return iter->second;
76  }
77  }
78  return 0;
79 }
80 
81 //_______________________________________________________________
83 {
84  unique_ptr<G4GDMLReadStructure> reader(new G4GDMLReadStructure());
85  G4GDMLParser gdmlParser(reader.get());
86  gdmlParser.SetOverlapCheck(OverlapCheck());
87  gdmlParser.Read(m_GDMPath, false);
88 
89  // alright the reader just puts everything into G4 Stores - endless fun
90  // print the show out, first solids:
91  // for (auto i=G4SolidStore::GetInstance()->begin(); i!=G4SolidStore::GetInstance()->end(); ++i)
92  // cout << "solid vol name: " << (*i)->GetName() << endl;
93  // for (auto i=G4LogicalVolumeStore::GetInstance()->begin(); i!=G4LogicalVolumeStore::GetInstance()->end(); i++)
94  // cout << "logvol name " << (*i)->GetName() << endl;
95 
97  for (set<string>::const_iterator its = mysubsys->assembly_iters().first; its != mysubsys->assembly_iters().second; ++its)
98  {
99  G4AssemblyVolume *avol = reader->GetAssembly(*its);
100  if (!avol)
101  {
102  cout << *its << "not found" << endl;
103  continue;
104  }
105  G4RotationMatrix *rotm = new G4RotationMatrix();
106  rotm->rotateX(m_Params->get_double_param("rot_x"));
107  rotm->rotateY(m_Params->get_double_param("rot_y"));
108  rotm->rotateZ(m_Params->get_double_param("rot_z"));
109  G4ThreeVector g4vec(m_Params->get_double_param("place_x"),
110  m_Params->get_double_param("place_y"),
111  m_Params->get_double_param("place_z"));
112  avol->MakeImprint(logicWorld, g4vec, rotm, 0, OverlapCheck());
113  delete rotm;
114  vector<G4VPhysicalVolume *>::iterator it = avol->GetVolumesIterator();
115  for (unsigned int i = 0; i < avol->TotalImprintedVolumes(); i++)
116  {
118  ++it;
119  }
120  }
121  for (set<string>::const_iterator its = mysubsys->logvol_iters().first; its != mysubsys->logvol_iters().second; ++its)
122  {
123  G4LogicalVolume *vol = reader->GetVolume(*its);
124  if (!vol)
125  {
126  cout << *its << "not found" << endl;
127  continue;
128  }
129 
130  G4RotationMatrix *rotm = new G4RotationMatrix();
131  rotm->rotateX(m_Params->get_double_param("rot_x"));
132  rotm->rotateX(m_Params->get_double_param("rot_y"));
133  rotm->rotateX(m_Params->get_double_param("rot_z"));
134  G4ThreeVector g4vec(m_Params->get_double_param("place_x"),
135  m_Params->get_double_param("place_y"),
136  m_Params->get_double_param("place_z"));
137  G4VPhysicalVolume *phys = new G4PVPlacement(rotm, g4vec,
138  vol,
139  G4String(GetName().c_str()),
140  logicWorld, false, 0, OverlapCheck());
142  }
143  AddHitNodes(topNode());
144  return;
145 }
146 
148 {
149  static int detid = -9999;
150  if (flag == insertassemblies)
151  {
152  // G4AssemblyVolumes naming convention:
153  // av_WWW_impr_XXX_YYY_ZZZ
154  // where:
155 
156  // WWW - assembly volume instance number
157  // XXX - assembly volume imprint number
158  // YYY - the name of the placed logical volume
159  // ZZZ - the logical volume index inside the assembly volume
160  // here we enter a new assembly volume and have to set the detector id which
161  // stays valid until we hit the next assembly volume
162  // the detector id is coded into YYY
163  if (physvol->GetName().find("av_") != string::npos && physvol->GetName().find("_impr_") != string::npos)
164  {
165  detid = -9999; // reset detid so we see if this is not handled here
166  std::vector<std::string> splitname;
167  boost::algorithm::split(splitname, physvol->GetName(), boost::is_any_of("_"));
168  if (splitname[4].find("AluStrips") != string::npos)
169  {
170  detid = 100;
171  }
172  else
173  {
174  string detprefix[] = {"VstStave", "FstContainerVolume", "BstContainerVolume", "Beampipe"};
175  //start layer count at 10 for VstStave
176  // at 20 for Fst
177  // at 30 for Bst
178  // at 40 for beampipe
179  int increase = 10;
180  for (auto toerase : detprefix)
181  {
182  size_t pos = splitname[4].find(toerase);
183  if (pos != string::npos)
184  {
185  detid = boost::lexical_cast<int>(splitname[4].erase(pos, toerase.length())) + increase;
186  break;
187  }
188  increase += 10;
189  }
190  }
191  }
192  if (detid < 0)
193  {
194  cout << "AllSiliconTrackerDetector::InsertVolumes detid is " << detid << " vol name " << physvol->GetName() << endl;
195  gSystem->Exit(1);
196  }
197  }
198  else
199  {
200  detid = 1;
201  }
202  G4LogicalVolume *logvol = physvol->GetLogicalVolume();
204  // cout << "mat: " << logvol->GetMaterial()->GetName() << endl;
205  // cout << "Adding " << physvol->GetName() << endl;
206  if (physvol->GetName().find("MimosaCore") != string::npos)
207  {
208  m_ActivePhysicalVolumesSet.insert(make_pair(physvol, detid));
209  }
210  else
211  {
212  if (m_AbsorberActive)
213  {
214  m_PassivePhysicalVolumesSet.insert(make_pair(physvol, -detid));
215  }
216  }
217  // G4 10.06 returns unsigned int for GetNoDaughters()
218  // lower version return int, need to cast to avoid compiler error
219  for (int i = 0; i < (int) logvol->GetNoDaughters(); ++i)
220  {
221  G4VPhysicalVolume *physvol = logvol->GetDaughter(i);
222  // here we decide which volumes are active
223  InsertVolumes(physvol, flag);
224  }
225  return;
226 }
227 
228 //_______________________________________________________________
229 void AllSiliconTrackerDetector::Print(const std::string &what) const
230 {
231  cout << "AllSiliconTracker Detector:" << endl;
232  if (what == "ALL" || what == "VOLUME")
233  {
234  cout << "Version 0.1" << endl;
235  cout << "Parameters:" << endl;
236  m_Params->Print();
237  }
238  return;
239 }
240 
241 int AllSiliconTrackerDetector::get_detid(const G4VPhysicalVolume *physvol, const int whichactive)
242 {
243  if (whichactive > 0)
244  {
245  return m_ActivePhysicalVolumesSet.find(physvol)->second;
246  }
247  else if (whichactive < 0)
248  {
249  return m_PassivePhysicalVolumesSet.find(physvol)->second;
250  }
251  else
252  {
253  cout << "AllSiliconTrackerDetector::get_detid invalid whichactive flag = 0" << endl;
254  gSystem->Exit(1);
255  }
256  return 0;
257 }
258 
260 {
261  PHNodeIterator iter(topNode);
262  PHCompositeNode *dstNode = dynamic_cast<PHCompositeNode *>(iter.findFirst("PHCompositeNode", "DST"));
263  PHNodeIterator dstIter(dstNode);
264  if (m_Params->get_int_param("active"))
265  {
266  map<string, int> nodes;
267  string myname;
268  if (SuperDetector() != "NONE")
269  {
270  myname = SuperDetector();
271  }
272  else
273  {
274  myname = GetName();
275  }
276  PHCompositeNode *DetNode = dynamic_cast<PHCompositeNode *>(dstIter.findFirst("PHCompositeNode", myname));
277  if (!DetNode)
278  {
279  DetNode = new PHCompositeNode(myname);
280  dstNode->addNode(DetNode);
281  }
282  set<int> DetidSet;
283  for (auto iter = m_ActivePhysicalVolumesSet.begin(); iter != m_ActivePhysicalVolumesSet.end(); ++iter)
284  {
285  DetidSet.insert(iter->second);
286  }
287  ostringstream g4hitnodeactive;
288  for (auto iter = DetidSet.begin(); iter != DetidSet.end(); ++iter)
289  {
290  g4hitnodeactive.str("");
291  string region;
292  if (*iter < 20)
293  {
294  region = "_CENTRAL_";
295  }
296  else if (*iter < 30)
297  {
298  region = "_FORWARD_";
299  }
300  else
301  {
302  region = "_BACKWARD_";
303  }
304  g4hitnodeactive << "G4HIT_" << myname << region << *iter;
305  nodes.insert(make_pair(g4hitnodeactive.str(), *iter));
306  }
307  if (m_Params->get_int_param("absorberactive"))
308  {
309  string g4hitnodename = "G4HIT_ABSORBER_" + myname;
310  nodes.insert(make_pair(g4hitnodename, -1));
311  }
312  for (auto nodename : nodes)
313  {
314  PHG4HitContainer *g4_hits = findNode::getClass<PHG4HitContainer>(DetNode, nodename.first);
315  if (!g4_hits)
316  {
317  g4_hits = new PHG4HitContainer(nodename.first);
318  DetNode->addNode(new PHIODataNode<PHObject>(g4_hits, nodename.first, "PHObject"));
319  }
320  m_HitContainerMap.insert(make_pair(nodename.second, g4_hits));
321  }
322  }
323  return;
324 }
325 
327 {
328  map<int, PHG4HitContainer *>::const_iterator iter = m_HitContainerMap.find(i);
329  if (iter != m_HitContainerMap.end())
330  {
331  return iter->second;
332  }
333  cout << "bad index: " << i << endl;
334  gSystem->Exit(1);
335  exit(1);
336 }