Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4_FGEM_EIC.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4_FGEM_EIC.C
1 
9 using namespace std;
10 
11 void
12 FGEM_Init()
13 {
14 
15 }
16 
17 void
18 FGEMSetup(PHG4Reco* g4Reco, const int N_Sector = 8, //
19  const double min_eta = 1.45 //
20  )
21 {
22 
23  const double tilt = .1;
24 
25  string name;
26  double etamax;
27  double etamin;
28  double zpos;
29  PHG4SectorSubsystem *gem;
30 
31  make_GEM_station("FGEM_0", g4Reco, 17, 1.01, 2.7, N_Sector);
32  make_GEM_station("FGEM_1", g4Reco, 62, 2.15, 4.0, N_Sector);
33 
35 
36  name = "FGEM_2";
37  etamax = 4;
38  etamin = min_eta;
39  zpos = 1.2e2;
40 
41  gem = new PHG4SectorSubsystem(name.c_str());
42 
43  gem->get_geometry().set_normal_polar_angle(tilt);
44  gem->get_geometry().set_normal_start(
45  zpos * PHG4Sector::Sector_Geometry::Unit_cm(), 0);
46  gem->get_geometry().set_min_polar_angle(
47  PHG4Sector::Sector_Geometry::eta_to_polar_angle(etamax));
48  gem->get_geometry().set_max_polar_angle(
49  PHG4Sector::Sector_Geometry::eta_to_polar_angle(etamin));
50  gem->get_geometry().set_max_polar_edge(
51  PHG4Sector::Sector_Geometry::FlatEdge());
52  gem->get_geometry().set_material("G4_METHANE");
53  gem->get_geometry().set_N_Sector(N_Sector);
54  gem->OverlapCheck(overlapcheck);
55  AddLayers_MiniTPCDrift(gem);
56  gem->get_geometry().AddLayers_HBD_GEM();
57  g4Reco->registerSubsystem(gem);
58 
60 
61  name = "FGEM_3";
62  etamax = 4;
63  etamin = min_eta;
64  zpos = 1.6e2;
65  gem = new PHG4SectorSubsystem(name.c_str());
66 
67  gem->SuperDetector(name);
68  gem->get_geometry().set_normal_polar_angle(tilt);
69  gem->get_geometry().set_normal_start(
70  zpos * PHG4Sector::Sector_Geometry::Unit_cm(), 0);
71  gem->get_geometry().set_min_polar_angle(
72  PHG4Sector::Sector_Geometry::eta_to_polar_angle(etamax));
73  gem->get_geometry().set_max_polar_angle(
74  PHG4Sector::Sector_Geometry::eta_to_polar_angle(2));
75  gem->get_geometry().set_max_polar_edge(
76  PHG4Sector::Sector_Geometry::FlatEdge());
77  gem->get_geometry().set_material("G4_METHANE");
78  gem->get_geometry().set_N_Sector(N_Sector);
79  gem->OverlapCheck(overlapcheck);
80  AddLayers_MiniTPCDrift(gem);
81  gem->get_geometry().AddLayers_HBD_GEM();
82  g4Reco->registerSubsystem(gem);
83 
84  gem = new PHG4SectorSubsystem(name + "_LowerEta");
85  gem->SuperDetector(name);
86 
87  zpos = zpos
88  - (zpos * sin(tilt)
89  + zpos * cos(tilt)
90  * tan(PHG4Sector::Sector_Geometry::eta_to_polar_angle(2) - tilt))
91  * sin(tilt);
92 
93  gem->get_geometry().set_normal_polar_angle(
94  (PHG4Sector::Sector_Geometry::eta_to_polar_angle(min_eta)
95  + PHG4Sector::Sector_Geometry::eta_to_polar_angle(2)) / 2);
96  gem->get_geometry().set_normal_start(
97  zpos * PHG4Sector::Sector_Geometry::Unit_cm(),
98  PHG4Sector::Sector_Geometry::eta_to_polar_angle(2));
99  gem->get_geometry().set_min_polar_angle(
100  PHG4Sector::Sector_Geometry::eta_to_polar_angle(2));
101  gem->get_geometry().set_max_polar_angle(
102  PHG4Sector::Sector_Geometry::eta_to_polar_angle(min_eta));
103  gem->get_geometry().set_material("G4_METHANE");
104  gem->get_geometry().set_N_Sector(N_Sector);
105  gem->get_geometry().set_min_polar_edge(
106  PHG4Sector::Sector_Geometry::FlatEdge());
107 
108  AddLayers_MiniTPCDrift(gem);
109  gem->get_geometry().AddLayers_HBD_GEM();
110  gem->OverlapCheck(overlapcheck);
111  g4Reco->registerSubsystem(gem);
112 
114 
115  name = "FGEM_4";
116  etamax = 4;
117  etamin = min_eta;
118  zpos = 2.75e2;
119  gem = new PHG4SectorSubsystem(name.c_str());
120 
121  gem->SuperDetector(name);
122  gem->get_geometry().set_normal_polar_angle(tilt);
123  gem->get_geometry().set_normal_start(
124  zpos * PHG4Sector::Sector_Geometry::Unit_cm(), 0);
125  gem->get_geometry().set_min_polar_angle(
126  PHG4Sector::Sector_Geometry::eta_to_polar_angle(etamax));
127  gem->get_geometry().set_max_polar_angle(
128  PHG4Sector::Sector_Geometry::eta_to_polar_angle(2));
129  gem->get_geometry().set_max_polar_edge(
130  PHG4Sector::Sector_Geometry::FlatEdge());
131  gem->get_geometry().set_material("G4_METHANE");
132  gem->get_geometry().set_N_Sector(N_Sector);
133  gem->OverlapCheck(overlapcheck);
134  AddLayers_MiniTPCDrift(gem);
135  gem->get_geometry().AddLayers_HBD_GEM();
136  g4Reco->registerSubsystem(gem);
137 
138  zpos = zpos
139  - (zpos * sin(tilt)
140  + zpos * cos(tilt)
141  * tan(PHG4Sector::Sector_Geometry::eta_to_polar_angle(2) - tilt))
142  * sin(tilt);
143 
144  gem = new PHG4SectorSubsystem(name + "_LowerEta");
145  gem->SuperDetector(name);
146 
147  gem->get_geometry().set_normal_polar_angle(
148  (PHG4Sector::Sector_Geometry::eta_to_polar_angle(min_eta)
149  + PHG4Sector::Sector_Geometry::eta_to_polar_angle(2)) / 2);
150  gem->get_geometry().set_normal_start(
151  zpos * PHG4Sector::Sector_Geometry::Unit_cm(),
152  PHG4Sector::Sector_Geometry::eta_to_polar_angle(2));
153  gem->get_geometry().set_min_polar_angle(
154  PHG4Sector::Sector_Geometry::eta_to_polar_angle(2));
155  gem->get_geometry().set_max_polar_angle(
156  PHG4Sector::Sector_Geometry::eta_to_polar_angle(min_eta));
157  gem->get_geometry().set_material("G4_METHANE");
158  gem->get_geometry().set_N_Sector(N_Sector);
159  gem->get_geometry().set_min_polar_edge(
160  PHG4Sector::Sector_Geometry::FlatEdge());
161 
162  AddLayers_MiniTPCDrift(gem);
163  gem->get_geometry().AddLayers_HBD_GEM();
164  gem->OverlapCheck(overlapcheck);
165  g4Reco->registerSubsystem(gem);
166 
168 
169 }
170 
172 void
173 AddLayers_MiniTPCDrift(PHG4SectorSubsystem *gem)
174 {
175  assert(gem);
176 
177  const double cm = PHG4Sector::Sector_Geometry::Unit_cm();
178  const double mm = .1 * cm;
179  const double um = 1e-3 * mm;
180 
181  // const int N_Layers = 70; // used for mini-drift TPC timing digitalization
182  const int N_Layers = 1; // simplified setup
183  const double thickness = 2 * cm;
184 
185  gem->get_geometry().AddLayer("EntranceWindow", "G4_MYLAR", 25 * um, false,
186  100);
187  gem->get_geometry().AddLayer("Cathode", "G4_GRAPHITE", 10 * um, false, 100);
188 
189  for (int d = 1; d <= N_Layers; d++)
190  {
191  stringstream s;
192  s << "DriftLayer_";
193  s << d;
194 
195  gem->get_geometry().AddLayer(s.str(), "G4_METHANE", thickness / N_Layers,
196  true);
197 
198  }
199 }
200 
201 int
202 make_GEM_station(string name, PHG4Reco* g4Reco, double zpos, double etamin,
203  double etamax, const int N_Sector = 8)
204 {
205 
206  // cout
207  // << "make_GEM_station - GEM construction with PHG4SectorSubsystem - make_GEM_station_EdgeReadout of "
208  // << name << endl;
209 
210  double polar_angle = 0;
211 
212  if (zpos < 0)
213  {
214  zpos = -zpos;
215  polar_angle = TMath::Pi();
216 
217  }
218  if (etamax < etamin)
219  {
220  double t = etamax;
221  etamax = etamin;
222  etamin = t;
223  }
224 
225  PHG4SectorSubsystem *gem;
226  gem = new PHG4SectorSubsystem(name.c_str());
227 
228  gem->SuperDetector(name);
229 
230  gem->get_geometry().set_normal_polar_angle(polar_angle);
231  gem->get_geometry().set_normal_start(
232  zpos * PHG4Sector::Sector_Geometry::Unit_cm());
233  gem->get_geometry().set_min_polar_angle(
234  PHG4Sector::Sector_Geometry::eta_to_polar_angle(etamax));
235  gem->get_geometry().set_max_polar_angle(
236  PHG4Sector::Sector_Geometry::eta_to_polar_angle(etamin));
237  gem->get_geometry().set_max_polar_edge(
238  PHG4Sector::Sector_Geometry::FlatEdge());
239  gem->get_geometry().set_min_polar_edge(
240  PHG4Sector::Sector_Geometry::FlatEdge());
241  gem->get_geometry().set_N_Sector(N_Sector);
242  gem->get_geometry().set_material("G4_METHANE");
243  gem->OverlapCheck(overlapcheck);
244 
245  AddLayers_MiniTPCDrift(gem);
246  gem->get_geometry().AddLayers_HBD_GEM();
247  g4Reco->registerSubsystem(gem);
248 
249 }
250 
251 void FGEM_FastSim_Reco(int verbosity = 0) {
252 
253  //---------------
254  // Load libraries
255  //---------------
256 
257  gSystem->Load("libfun4all.so");
258  gSystem->Load("libg4hough.so");
259 
260  //---------------
261  // Fun4All server
262  //---------------
263 
264  Fun4AllServer *se = Fun4AllServer::instance();
265 
266  PHG4TrackFastSim* kalman = new PHG4TrackFastSim("PHG4TrackFastSim");
267  kalman->Verbosity(0);
268 
269  kalman->set_use_vertex_in_fitting(true);
270  kalman->set_vertex_xy_resolution(50E-4);
271  kalman->set_vertex_z_resolution(50E-4);
272 
273  kalman->set_detector_type(PHG4TrackFastSim::Vertical_Plane); // Vertical_Plane, Cylinder
274  kalman->set_phi_resolution(50E-4);
275  kalman->set_r_resolution(1.);
276 
277  kalman->set_pat_rec_hit_finding_eff(1.);
278  kalman->set_pat_rec_noise_prob(0.);
279 
280  std::string phg4hits_names[] = {"G4HIT_FGEM_0","G4HIT_FGEM_1","G4HIT_FGEM_2","G4HIT_FGEM_3","G4HIT_FGEM_4"};
281  kalman->set_phg4hits_names(phg4hits_names, 5);
282  kalman->set_sub_top_node_name("SVTX");
283  kalman->set_trackmap_out_name("SvtxTrackMap_FastSimEtaPlus");
284 
285  // Saved track states (projections)
286  std::string state_names[] = {"FEMC","FHCAL"};
287  kalman->set_state_names(state_names, 2);
288 
289  kalman->set_fit_alg_name("KalmanFitterRefTrack");//
290  kalman->set_primary_assumption_pid(13);
291  kalman->set_do_evt_display(false);
292 
293  se->registerSubsystem(kalman);
294 
295 }