Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
SetupDualRICHAnalyzer.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file SetupDualRICHAnalyzer.cc
1 // Some documentation, possibly //
2 
3 #include "TMath.h"
4 
6 #include "dualrich_analyzer.h"
7 
8 using namespace std;
9 
10 
12 {
13  _analyzer = new eic_dual_rich();
14 }
15 
16 double
17 SetupDualRICHAnalyzer::calculate_emission_angle( double m_emi[3], double momv[3], PHG4Hit * hit_i )
18 {
19  /* Input parameters for indirect ray tracing algorithm */
20  double Ex = m_emi[0];
21  double Ey = m_emi[1];
22  double Ez = m_emi[2];
23 
24  double Dx = hit_i->get_x(0);
25  double Dy = hit_i->get_y(0);
26  double Dz = hit_i->get_z(0);
27 
28  double vx = momv[0];
29  double vy = momv[1];
30  double vz = momv[2];
31 
32  int sec = hit_i->get_detid();
33  double cx = -18.5*TMath::Sin(sec*TMath::Pi()/4); // mirror center of each octant
34  double cy = 18.5*TMath::Cos(sec*TMath::Pi()/4);
35  double cz = 75;
36 
37  int select_radiator=0;
38 
39  /* Set mirror parameters */
40  double R_mirror = 195; // cm
41  _analyzer->set_mirror(cx, cy, cz, R_mirror);
42 
43  /* Call algorithm to determine emission angle of photon i w.r.t. track j */
44  float theta_c = _analyzer->ind_ray(Ex, Ey, Ez, Dx, Dy, Dz, vx, vy, vz, select_radiator); //Indirect Ray Tracing
45 
46  return theta_c;
47 }
48 
49 bool
51 {
52  arr_mom[0] = 0;
53  arr_mom[1] = 0;
54  arr_mom[2] = 0;
55 
56  /* Get the track's particle, then get truth momentum */
57  PHG4Particle* particle = truthinfo->GetParticle( track->get_truth_track_id() );
58  if (particle){
59  arr_mom[0] = particle->get_px();
60  arr_mom[1] = particle->get_py();
61  arr_mom[2] = particle->get_pz();
62  return true;
63  }
64 
65  return false;
66 }
67 
68 bool
70 {
71  arr_mom[0] = 0;
72  arr_mom[1] = 0;
73  arr_mom[2] = 0;
74 
75  /* Get photon emission points, then get truth momentum using the vector from first to last emitted photon */
76  vector<double> emix;
77  vector<double> emiy;
78  vector<double> emiz;
79 
80  /* Loop over all G4Hits in container (i.e. RICH photons in event) */
81  PHG4HitContainer::ConstRange rich_hits_begin_end = richhits->getHits();
82  PHG4HitContainer::ConstIterator rich_hits_iter;
83 
84  for (rich_hits_iter = rich_hits_begin_end.first; rich_hits_iter != rich_hits_begin_end.second; ++rich_hits_iter)
85  {
86  PHG4Hit *hit_i = rich_hits_iter->second;
87  PHG4Particle *particle = truthinfo->GetParticle( hit_i->get_trkid() );
88  PHG4VtxPoint *vertex = truthinfo->GetVtx( particle->get_vtx_id() );
89 
90  /* Should eventually compare that vertices are along the correct tracks, for now we run single tracks */
91  emix.push_back( vertex->get_x() );
92  emiy.push_back( vertex->get_y() );
93  emiz.push_back( vertex->get_z() );
94  }
95 
96 
97  vector<double>::iterator first;
98  vector<double>::iterator last;
99  double dx=0;
100  double dy=0;
101  double dz=0;
102 
103  /* Use first-to-last (or first-to-11th) */
104  first = std::min_element(emiz.begin(),emiz.end());
105  //last = std::max_element(emiz.begin(),emiz.end());
106  double p1 = std::distance(emiz.begin(),first);
107  //double p2 = std::distance(emiz.begin(),last);
108  if (emiz.size() > p1+11){
109  dx = emix.at(p1+11) - emix.at(p1);
110  dy = emiy.at(p1+11) - emiy.at(p1);
111  dz = emiz.at(p1+11) - emiz.at(p1);
112  }
113 
114  /* Fill "momentum" */
115  arr_mom[0] = dx;
116  arr_mom[1] = dy;
117  arr_mom[2] = dz;
118 
119  return true;
120 }