Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Plot-Measured_Energy_Over_True_Energy-EMC.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file Plot-Measured_Energy_Over_True_Energy-EMC.C
1 #include <iostream>
2 #include <cstdlib>
3 #include <memory>
4 #include <sstream>
5 #include <string>
6 #include <cstring>
7 
8 #include "TFile.h"
9 #include "TTree.h"
10 #include "TCanvas.h"
11 #include "TLegend.h"
12 
13 #define NELEMS(arr) (sizeof(arr)/sizeof(arr[0]))
14 
15 /* Directory where data is stored for plots */
16 const char *const data_directory =
17  "/sphenix/user/gregtom3/data/Summer2018/ECAL_energy_studies";
18 /* The energy levels we have in GeV */
19 const static int energy_levels[] = { 1, 2, 5, 10, 20 };
20 
21 /* The particle types we have */
22 enum particle_type { electron, pion };
23 /* The detectors we have */
24 enum detector { cemc, eemc, femc };
25 
26 TTree *load_tree(const char *const file_name, const char *const tree_name);
27 void fill_histogram(TH1F * const h, TTree * const t, const Float_t true_energy,
28  const Float_t min_value, const bool normalize);
29 void histogram_to_png(TH1F * const h_pion_CEMC, TH1F * const h_electron_CEMC,
30  TH1F * const h_pion_EEMC, TH1F * const h_electron_EEMC,
31  TH1F * const h_pion_FEMC, TH1F * const h_electron_FEMC,
32  const char *const title, const char *const save_file_name,
33  const char *const cemc_label,
34  const char *const eemc_label,
35  const char *const femc_label);
36 char *generate_name(const particle_type p, const int particle_energy_gev,
37  const detector d);
38 char *generate_file_path(const particle_type p, const int particle_energy_gev,
39  const detector d);
40 char *generate_save_name(const int particle_energy_gev);
41 char *generate_title(const int particle_energy_gev);
42 char *generate_label(const int particle_energy_gev, const detector d);
43 
44 void Plot_Measured_Energy_Over_True_Energy_EMC()
45 {
46 
47  /*
48  * sPHENIX Style
49  */
50 
51  gROOT->LoadMacro
52  ("/sphenix/user/gregtom3/SBU/research/macros/macros/sPHENIXStyle/sPhenixStyle.C");
53  SetsPhenixStyle();
54 
55  gROOT->SetBatch(kTRUE);
56 
57  /*
58  * Base Histogram (Recreated from Matching Plots)
59  */
60 
61  TH1F *h_base = new TH1F("h_base", "", 25, 0.0, 2.5);
62  TH1F *h_base_e = (TH1F *) h_base->Clone();
63  TH1F *h_base_p = (TH1F *) h_base->Clone();
64  h_base_e->SetLineColor(kRed);
65  h_base_p->SetLineColor(kBlue);
66 
67  /* iterate through energy levels and create plots for Pions and Electrons
68  * in the CEMC and EEMC detectors */
69  for (size_t i = 0; i < NELEMS(energy_levels); ++i) {
70  /* CEMC */
71  TH1F *const h_pion_cemc = h_base_p->Clone();
72  TH1F *const h_electron_cemc = h_base_e->Clone();
73  h_pion_cemc->SetName(generate_name
74  (pion, energy_levels[i], cemc));
75  h_electron_cemc->SetName(generate_name
76  (electron, energy_levels[i], cemc));
77 
78  TTree *const t_pion_cemc =
79  load_tree(generate_file_path(pion, energy_levels[i], cemc),
80  "ntp_cluster");
81  TTree *const t_electron_cemc =
82  load_tree(generate_file_path
83  (electron, energy_levels[i], cemc),
84  "ntp_cluster");
85 
86  fill_histogram(h_pion_cemc, t_pion_cemc, energy_levels[i], 0.3,
87  true);
88  fill_histogram(h_electron_cemc, t_electron_cemc,
89  energy_levels[i], 0.3, true);
90 
91  /* EEMC */
92  TH1F *const h_pion_eemc = h_base_p->Clone();
93  TH1F *const h_electron_eemc = h_base_e->Clone();
94  h_pion_eemc->SetName(generate_name
95  (pion, energy_levels[i], eemc));
96  h_electron_eemc->SetName(generate_name
97  (electron, energy_levels[i], eemc));
98 
99  TTree *const t_pion_eemc =
100  load_tree(generate_file_path(pion, energy_levels[i], eemc),
101  "ntp_cluster");
102  TTree *const t_electron_eemc =
103  load_tree(generate_file_path
104  (electron, energy_levels[i], eemc),
105  "ntp_cluster");
106 
107  fill_histogram(h_pion_eemc, t_pion_eemc, energy_levels[i], 0.3,
108  true);
109  fill_histogram(h_electron_eemc, t_electron_eemc,
110  energy_levels[i], 0.3, true);
111 
112  /* FEMC */
113  TH1F *const h_pion_femc = h_base_p->Clone();
114  TH1F *const h_electron_femc = h_base_e->Clone();
115  h_pion_femc->SetName(generate_name
116  (pion, energy_levels[i], femc));
117  h_electron_femc->SetName(generate_name
118  (electron, energy_levels[i], femc));
119 
120  TTree *const t_pion_femc =
121  load_tree(generate_file_path(pion, energy_levels[i], femc),
122  "ntp_cluster");
123  TTree *const t_electron_femc =
124  load_tree(generate_file_path
125  (electron, energy_levels[i], femc),
126  "ntp_cluster");
127 
128  fill_histogram(h_pion_femc, t_pion_femc, energy_levels[i], 0.3,
129  true);
130  fill_histogram(h_electron_femc, t_electron_femc,
131  energy_levels[i], 0.3, true);
132 
133  histogram_to_png(h_pion_cemc, h_electron_cemc,
134  h_pion_eemc, h_electron_eemc,
135  h_pion_femc, h_electron_femc,
136  generate_title(energy_levels[i]),
137  generate_save_name(energy_levels[i]),
138  generate_label(energy_levels[i], cemc),
139  generate_label(energy_levels[i], eemc),
140  generate_label(energy_levels[i], femc));
141  }
142 
143 }
144 
145 TTree *load_tree(const char *const file_name, const char *const tree_name)
146 {
147  return (TTree *) (new TFile(file_name, "READ"))->Get(tree_name);
148 }
149 
150 /*
151  * This function, when fed a histogram, tree,
152  * a floating point min_value variable, and a boolean, will fill each entry
153  * inside the specificed branch into the histogram, so long as each entry
154  * is a floating point number GREATER than the min_value. If normalize is
155  * set to 'true', the histogram will be normalize based on its number of
156  * entries. The axes titles are furthermore assumed to be generic and have
157  * been already set.
158  */
159 void fill_histogram(TH1F * const h, TTree * const t, const Float_t true_energy,
160  const Float_t min_value, const bool normalize)
161 {
162  Float_t measured_energy;
163 // Float_t true_energy;
164  t->SetBranchAddress("e", &measured_energy);
165 // t->SetBranchAddress("ge", &true_energy);
166  Float_t true_eta;
167  t->SetBranchAddress("geta", &true_eta);
168 
169  Int_t nentries = Int_t(t->GetEntries());
170 
171  for (Int_t i = 0; i < nentries; ++i) {
172  if (t->LoadTree(i) < 0)
173  break;
174 
175  t->GetEntry(i);
176  if (((true_eta > -0.5 && true_eta < 0.5)
177  || (true_eta > -3 && true_eta < -2)
178  || (true_eta > 2 && true_eta < 3))
179  && (measured_energy > min_value && true_energy > 0.1))
180  h->Fill(measured_energy / true_energy);
181  }
182  if (normalize)
183  h->Scale(1 / h->GetEntries());
184 
185  h->SetXTitle("E_{cluster} / E_{true}");
186  h->SetYTitle("entries / #scale[0.5]{#sum} entries ");
187 }
188 
189 void histogram_to_png(TH1F * const h_pion_CEMC, TH1F * const h_electron_CEMC,
190  TH1F * const h_pion_EEMC, TH1F * const h_electron_EEMC,
191  TH1F * const h_pion_FEMC, TH1F * const h_electron_FEMC,
192  const char *const title, const char *const save_file_name,
193  const char *const cemc_label,
194  const char *const eemc_label,
195  const char *const femc_label)
196 {
197  /* 3 for CEMC, EEMC, and FEMC */
198  const unsigned ndetectors = 3;
199  TCanvas cPNG("cPNG", title, gStyle->GetCanvasDefW() * 3,
200  gStyle->GetCanvasDefH());
201  TImage *img = TImage::Create();
202 
203  cPNG.Divide(3, 1);
204  cPNG.cd(1);
205  h_pion_CEMC->GetYaxis()->SetRangeUser(0.0001, 1);
206  h_electron_CEMC->GetYaxis()->SetRangeUser(0.0001, 1);
207  gPad->SetLogy();
208  h_pion_CEMC->Draw();
209  h_electron_CEMC->Draw("SAME");
210  gPad->RedrawAxis();
211 
212  auto cemc_legend = new TLegend(0.70, 0.9, 0.95, 0.65, cemc_label);
213  cemc_legend->AddEntry(h_pion_CEMC, "Pions", "l");
214  cemc_legend->AddEntry(h_electron_CEMC, "Electrons", "l");
215  cemc_legend->Draw();
216 
217  cPNG.cd(2);
218  h_pion_EEMC->GetYaxis()->SetRangeUser(0.0001, 1);
219  h_electron_EEMC->GetYaxis()->SetRangeUser(0.0001, 1);
220  gPad->SetLogy();
221  h_pion_EEMC->Draw();
222  h_electron_EEMC->Draw("SAME");
223  gPad->RedrawAxis();
224 
225  auto eemc_legend = new TLegend(0.65, 0.9, 0.9, 0.65, eemc_label);
226  eemc_legend->AddEntry(h_pion_EEMC, "Pions", "l");
227  eemc_legend->AddEntry(h_electron_EEMC, "Electrons", "l");
228  eemc_legend->Draw();
229 
230  cPNG.cd(3);
231  h_pion_FEMC->GetYaxis()->SetRangeUser(0.0001, 1);
232  h_electron_FEMC->GetYaxis()->SetRangeUser(0.0001, 1);
233  gPad->SetLogy();
234  h_pion_FEMC->Draw();
235  h_electron_FEMC->Draw("SAME");
236  gPad->RedrawAxis();
237 
238  auto femc_legend = new TLegend(0.65, 0.9, 0.9, 0.65, femc_label);
239  femc_legend->AddEntry(h_pion_FEMC, "Pions", "l");
240  femc_legend->AddEntry(h_electron_FEMC, "Electrons", "l");
241  femc_legend->Draw();
242 
243  img->FromPad(&cPNG);
244  img->WriteImage(save_file_name);
245 
246  delete img;
247 }
248 
249 char *strdup(const char *s)
250 {
251  char *const t = new char[strlen(s) + 1];
252  return strcpy(t, s);
253 }
254 
255 char *generate_name(const particle_type p, const int particle_energy_gev,
256  const detector d)
257 {
258  std::stringstream name;
259  name << "Histogram";
260  switch (p) {
261  case electron:
262  name << "_Electron";
263  break;
264  case pion:
265  name << "_Pion";
266  break;
267  }
268 
269  name << "_" << particle_energy_gev << "GeV";
270  switch (d) {
271  case cemc:
272  name << "_CEMC";
273  break;
274  case eemc:
275  name << "_EEMC";
276  break;
277  }
278 
279  return strdup(name.str().c_str());
280 }
281 
282 char *generate_file_path(const particle_type p,
283  const int particle_energy_gev, const detector d)
284 {
285  std::stringstream path;
286  path << data_directory;
287 
288  switch (p) {
289  case electron:
290  path << "/Electrons/Electrons";
291  break;
292  case pion:
293  path << "/Pions/Pions";
294  break;
295  }
296 
297  path << particle_energy_gev;
298  switch (d) {
299  case cemc:
300  path << "C";
301  break;
302  case eemc:
303  path << "E";
304  break;
305  case femc:
306  path << "F";
307  break;
308  }
309 
310  path << ".root";
311 
312  return strdup(path.str().c_str());
313 }
314 
315 char *generate_save_name(const int particle_energy_gev)
316 {
317  std::stringstream name;
318  name << "Electron-Pion-" << particle_energy_gev <<
319  "GeV-CEMC-EEMC-FEMC.png";
320 
321  return strdup(name.str().c_str());
322 }
323 
324 char *generate_title(const int particle_energy_gev)
325 {
326  std::stringstream title;
327  title << particle_energy_gev << "GeV";
328 
329  return strdup(title.str().c_str());
330 }
331 
332 char *generate_label(const int particle_energy_gev, const detector d)
333 {
334  std::stringstream label;
335 
336  switch (d) {
337  case cemc:
338  label << "CEMC ";
339  break;
340  case eemc:
341  label << "EEMC ";
342  break;
343  case femc:
344  label << "FEMC ";
345  break;
346  }
347  label << particle_energy_gev << "GeV";
348 
349  return strdup(label.str().c_str());
350 }