Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
TrackerEventDisplay.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file TrackerEventDisplay.cc
1 #include "TrackerEventDisplay.h"
2 
3 #include <trackbase/ActsGeometry.h>
4 #include <trackbase/TrkrCluster.h>
5 #include <trackbase/TrkrDefs.h>
6 #include <trackbase/TrkrHit.h>
7 #include <trackbase/TrkrHitSetContainer.h>
8 #include <trackbase/TrackFitUtils.h>
9 
10 #include <trackbase_historic/TrackSeedContainer_v1.h>
11 
12 #include <trackbase/TpcDefs.h>
13 
14 #include <trackbase/TrkrClusterContainer.h>
15 #include <trackbase/TrkrClusterHitAssoc.h>
16 #include <trackbase/TrkrClusterIterationMapv1.h>
17 #include <trackbase/TrkrHitSet.h>
18 
19 #include <g4detectors/PHG4TpcCylinderGeom.h>
20 #include <g4detectors/PHG4TpcCylinderGeomContainer.h>
21 
22 #include <fun4all/Fun4AllReturnCodes.h>
23 #include <fun4all/SubsysReco.h>
24 
25 #include <phool/getClass.h>
26 #include <phool/phool.h>
27 #include <phool/recoConsts.h>
28 
29 #include <TVector3.h>
30 #include <boost/format.hpp>
31 #include <boost/math/special_functions/sign.hpp>
32 
33 #include <cmath>
34 #include <iomanip>
35 #include <iostream>
36 #include <iterator>
37 #include <map>
38 #include <memory> // for shared_ptr
39 #include <set> // for _Rb_tree_cons...
40 #include <utility>
41 #include <vector>
42 
43 using namespace std;
44 
45 TrackerEventDisplay::TrackerEventDisplay(const string& /*name*/, const string& filename, const string& runnumber, const string& date)
46  : SubsysReco("TrackerEventDisplay")
47  , _hit(true)
48  , _cluster(false)
49  , _ievent(0)
50  , _filename(filename)
51  , _runnumber(runnumber)
52  , _date(date)
53 {
54 }
55 
56 TrackerEventDisplay::~TrackerEventDisplay()
57 {
58 }
59 
60 int TrackerEventDisplay::Init(PHCompositeNode* /*topNode*/)
61 {
62  _ievent = 0;
63 
64  return Fun4AllReturnCodes::EVENT_OK;
65 }
66 
67 int TrackerEventDisplay::InitRun(PHCompositeNode*)
68 {
69  return Fun4AllReturnCodes::EVENT_OK;
70 }
71 
72 int TrackerEventDisplay::process_event(PHCompositeNode* topNode)
73 {
74  makeJsonFile(topNode);
75  ++_ievent;
76  return Fun4AllReturnCodes::EVENT_OK;
77 }
78 
79 int TrackerEventDisplay::End(PHCompositeNode* /*topNode*/)
80 {
81  if (Verbosity() > 1)
82  {
83  cout << "========================= TrackerEventDisplay::End() ============================" << endl;
84  cout << " " << _ievent << " events of output written to: " << _filename << endl;
85  cout << "===========================================================================" << endl;
86  }
87 
88  return Fun4AllReturnCodes::EVENT_OK;
89 }
90 
91 void TrackerEventDisplay::makeJsonFile(PHCompositeNode* topNode)
92 {
93  if (Verbosity() > 1)
94  {
95  cout << "TrackerEventDisplay::makeJsonFile() entered" << endl;
96  }
97 
98  ActsGeometry* tgeometry = findNode::getClass<ActsGeometry>(topNode, "ActsGeometry");
99  if (!tgeometry)
100  {
101  std::cout << "No Acts geometry on node tree. Can't continue."
102  << std::endl;
103  return;
104  }
105 
106  PHG4TpcCylinderGeomContainer* geom_container =
107  findNode::getClass<PHG4TpcCylinderGeomContainer>(topNode, "CYLINDERCELLGEOM_SVTX");
108  if (!geom_container)
109  {
110  std::cout << PHWHERE << "ERROR: Can't find node CYLINDERCELLGEOM_SVTX" << std::endl;
111  return;
112  }
113 
114  //--------------------
115  // fill the Hit json
116  //--------------------
117 
118  if (_hit)
119  {
120  bool firstHit = true;
121 
122  outdata.open((_filename+"_event"+std::to_string(_ievent)+"_hits.json").c_str(),std::ofstream::out | std::ofstream::trunc);
123 
124  if( !outdata )
125  { // file couldn't be opened
126  cerr << "ERROR: file could not be opened" << endl;
127  exit(1);
128  }
129 
130  outdata << "{\n \"EVENT\": {\n \"runid\":" << _runnumber << ", \n \"evtid\": 1, \n \"time\": 0, \n \"timeStr\": \"2023-08-23, 15:23:30 EST\", \n \"type\": \"Cosmics\", \n \"s_nn\": 0, \n \"B\": 0.0,\n \"pv\": [0,0,0],\n \"runstats\": [ \n \"sPHENIX Time Projection Chamber\", \"" << _date << ", Run " << _runnumber << " - Event " << _ievent << "\", \"All Hits in Event\"] \n },\n" << endl;
131 
132  outdata << " \"META\": {\n \"HITS\": {\n \"INNERTRACKER\": {\n \"type\": \"3D\",\n \"options\": {\n \"size\": 2,\n \"color\": 16777215\n } \n },\n" << endl;
133  outdata << " \"TRACKHITS\": {\n \"type\": \"3D\",\n \"options\": {\n \"size\": 2,\n \"transparent\": 0.5,\n \"color\": 16777215\n } \n },\n" << endl;
134  outdata << " \"JETS\": {\n \"type\": \"JET\",\n \"options\": {\n \"rmin\": 0,\n \"rmax\": 78,\n \"emin\": 0,\n \"emax\": 30,\n \"color\": 16777215,\n \"transparent\": 0.5 \n }\n }\n }\n }\n," << endl;
135  outdata << " \"HITS\": {\n \"CEMC\":[{\"eta\": 0, \"phi\": 0, \"e\": 0}\n ],\n \"HCALIN\": [{\"eta\": 0, \"phi\": 0, \"e\": 0}\n ],\n \"HCALOUT\": [{\"eta\": 0, \"phi\": 0, \"e\": 0}\n \n ],\n\n" << endl;
136  outdata << " \"TRACKHITS\": [\n\n ";
137 
138  auto m_tGeometry = findNode::getClass<ActsGeometry>(topNode, "ActsGeometry");
139 
140  if (Verbosity() >= 1)
141  {
142  cout << "Filling hit json" << endl;
143  }
144  // need things off of the DST...
145  TrkrHitSetContainer* hitmap = findNode::getClass<TrkrHitSetContainer>(topNode, "TRKR_HITSET");
146 
147  if (hitmap)
148  {
149  TrkrHitSetContainer::ConstRange all_hitsets = hitmap->getHitSets();
150  for (TrkrHitSetContainer::ConstIterator iter = all_hitsets.first;
151  iter != all_hitsets.second;
152  ++iter)
153  {
154  const TrkrDefs::hitsetkey hitset_key = iter->first;
155  TrkrHitSet* hitset = iter->second;
156  // get all hits for this hitset
157  TrkrHitSet::ConstRange hitrangei = hitset->getHits();
158  for (TrkrHitSet::ConstIterator hitr = hitrangei.first;
159  hitr != hitrangei.second;
160  ++hitr)
161  {
162  TrkrDefs::hitkey hit_key = hitr->first;
163  TrkrHit* hit = hitr->second;
164  //float event = _ievent;
165  //float hitID = hit_key;
166  //float e = hit->getEnergy();
167  float adc = hit->getAdc();
168  float layer_local = TrkrDefs::getLayer(hitset_key);
169  //float sector = TpcDefs::getSectorId(hitset_key);
170  float side = TpcDefs::getSide(hitset_key);
171  //float cellID = 0;
172  //float ecell = hit->getAdc();
173 
174  float phibin = NAN;
175  float tbin = NAN;
176  //float phi = NAN;
177  float phi_center = NAN;
178  float x = NAN;
179  float y = NAN;
180  float z = NAN;
181 
182  if (TrkrDefs::getTrkrId(hitset_key) == TrkrDefs::TrkrId::tpcId)
183  {
184  PHG4TpcCylinderGeom* GeoLayer_local = geom_container->GetLayerCellGeom(layer_local);
185  double radius = GeoLayer_local->get_radius();
186  phibin = (float) TpcDefs::getPad(hit_key);
187  tbin = (float) TpcDefs::getTBin(hit_key);
188  //phi = GeoLayer_local->get_phicenter(phibin);
189 
190  double zdriftlength = tbin * m_tGeometry->get_drift_velocity() * AdcClockPeriod;
191  // convert z drift length to z position in the TPC
192  // cout << " tbin: " << tbin << " vdrift " <<m_tGeometry->get_drift_velocity() << " l drift: " << zdriftlength <<endl;
193  unsigned short NTBins = (unsigned short) GeoLayer_local->get_zbins();
194  double m_tdriftmax = AdcClockPeriod * NTBins / 2.0;
195  double clusz = (m_tdriftmax * m_tGeometry->get_drift_velocity()) - zdriftlength;
196  if (side == 0)
197  {
198  clusz = -clusz;
199  }
200  z = clusz;
201  phi_center = GeoLayer_local->get_phicenter(phibin);
202  x = radius * cos(phi_center);
203  y = radius * sin(phi_center);
204 
205  stringstream spts;
206 
207  if (firstHit) firstHit = false;
208  else spts << ",";
209 
210  spts << "{ \"x\": ";
211  spts << x;
212  spts << ", \"y\": ";
213  spts << y;
214  spts << ", \"z\": ";
215  spts << z;
216  spts << ", \"e\": ";
217  spts << adc;
218  spts << "}";
219 
220  outdata << (boost::format("%1%") % spts.str());
221  spts.clear();
222  spts.str("");
223  }
224  }
225  }
226  }
227  outdata << "],\n \"JETS\": [\n ]\n }," << endl;
228  outdata << "\"TRACKS\": {" << endl;
229  outdata <<"\""<<"INNERTRACKER"<<"\": [";
230  outdata << "]" << endl;
231  outdata << "}" << endl;
232  outdata << "}" << endl;
233  outdata.close();
234  }
235 
236  //------------------------
237  // fill the Cluster JSON
238  //------------------------
239 
240  if (_cluster)
241  {
242  bool firstHit = true;
243 
244  outdata.open((_filename+"_event"+std::to_string(_ievent)+"_clusters.json").c_str(),std::ofstream::out | std::ofstream::trunc);
245 
246  if( !outdata )
247  { // file couldn't be opened
248  cerr << "ERROR: file could not be opened" << endl;
249  exit(1);
250  }
251 
252  outdata << "{\n \"EVENT\": {\n \"runid\":" << _runnumber << ", \n \"evtid\": 1, \n \"time\": 0, \n \"timeStr\": \"2023-08-23, 15:23:30 EST\", \n \"type\": \"Cosmics\", \n \"s_nn\": 0, \n \"B\": 0.0,\n \"pv\": [0,0,0],\n \"runstats\": [ \n \"sPHENIX Time Projection Chamber\", \"" << _date << ", Run " << _runnumber << " - Event " << _ievent << "\", \"All Clusters in Event\"] \n },\n" << endl;
253 
254  outdata << " \"META\": {\n \"HITS\": {\n \"INNERTRACKER\": {\n \"type\": \"3D\",\n \"options\": {\n \"size\": 2,\n \"color\": 16777215\n } \n },\n" << endl;
255  outdata << " \"TRACKHITS\": {\n \"type\": \"3D\",\n \"options\": {\n \"size\": 2,\n \"transparent\": 0.5,\n \"color\": 16777215\n } \n },\n" << endl;
256  outdata << " \"JETS\": {\n \"type\": \"JET\",\n \"options\": {\n \"rmin\": 0,\n \"rmax\": 78,\n \"emin\": 0,\n \"emax\": 30,\n \"color\": 16777215,\n \"transparent\": 0.5 \n }\n }\n }\n }\n," << endl;
257  outdata << " \"HITS\": {\n \"CEMC\":[{\"eta\": 0, \"phi\": 0, \"e\": 0}\n ],\n \"HCALIN\": [{\"eta\": 0, \"phi\": 0, \"e\": 0}\n ],\n \"HCALOUT\": [{\"eta\": 0, \"phi\": 0, \"e\": 0}\n \n ],\n\n" << endl;
258  outdata << " \"TRACKHITS\": [\n\n ";
259 
260  if (Verbosity() > 1)
261  {
262  cout << "Filling cluster json " << endl;
263  }
264 
265  // need things off of the DST...
266  TrkrClusterContainer* clustermap = findNode::getClass<TrkrClusterContainer>(topNode, "CORRECTED_TRKR_CLUSTER");
267  if (!clustermap)
268  {
269  clustermap = findNode::getClass<TrkrClusterContainer>(topNode, "TRKR_CLUSTER");
270  }
271 
272  if (Verbosity() > 1)
273  {
274  if (clustermap != nullptr)
275  {
276  cout << "got clustermap" << endl;
277  }
278  else
279  {
280  cout << "no clustermap" << endl;
281  }
282  }
283 
284  if (clustermap)
285  {
286  for (const auto& hitsetkey : clustermap->getHitSetKeys())
287  {
288  auto range = clustermap->getClusters(hitsetkey);
289  for (auto iter = range.first; iter != range.second; ++iter)
290  {
291  TrkrDefs::cluskey cluster_key = iter->first;
292  TrkrCluster* cluster = clustermap->findCluster(cluster_key);
293 
294  Acts::Vector3 cglob;
295  cglob = tgeometry->getGlobalPosition(cluster_key, cluster);
296  float x = cglob(0);
297  float y = cglob(1);
298  float z = cglob(2);
299  float adc = cluster->getAdc();
300 
301  stringstream spts;
302 
303  if (firstHit) firstHit = false;
304  else spts << ",";
305 
306  spts << "{ \"x\": ";
307  spts << x;
308  spts << ", \"y\": ";
309  spts << y;
310  spts << ", \"z\": ";
311  spts << z;
312  spts << ", \"e\": ";
313  spts << adc;
314  spts << "}";
315 
316  outdata << (boost::format("%1%") % spts.str());
317  spts.clear();
318  spts.str("");
319  }
320  }
321  }
322  outdata << "],\n \"JETS\": [\n ]\n }," << endl;
323  outdata << "\"TRACKS\": {" << endl;
324  outdata <<"\""<<"INNERTRACKER"<<"\": [";
325  outdata << "]" << endl;
326  outdata << "}" << endl;
327  outdata << "}" << endl;
328  outdata.close();
329  }
330  return;
331 }