Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
plot_true_false_positive.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file plot_true_false_positive.C
1 #include "TGraph.h"
2 vector< TString > v_cuts_CEMC;
3 vector< TString > v_cuts_CEMC_short;
4 vector< TString > v_cuts_EEMC;
5 vector< TString > v_cuts_EEMC_short;
6 vector< TString > v_cuts_FEMC;
7 vector< TString > v_cuts_FEMC_short;
8 
9 void add_CEMC_Cut(double, double);
10 
11 void add_EEMC_Cut(double);
12 
13 void add_FEMC_Cut(double);
14 
16 {
17 
18  /*----------------------------------------- */
19  /* Add Cuts Here */
20  /* */
21  /* add_CEMC_cut(float_1, float_2) */
22  /* Cut 1 : E/p > float_1 */
23  /* Cut 2 : prob > float_2 */
24  /* (if float_#=0, Cut # is not applied */
25  /* add_EEMC_cut(float_1) */
26  /* Cut 1 : E/True E > float_1 */
27  /* add_FEMC_cut(float_1) */
28  /* Cut 1 : E/True E > float_1 */
29  /* ---------------------------------------- */
30 
31  add_CEMC_Cut(0.807,0.02);
32  add_CEMC_Cut(0.703,0.02);
33 
34  add_EEMC_Cut(0.703);
35  add_EEMC_Cut(0.807);
36  add_EEMC_Cut(0.95);
37 
38  add_FEMC_Cut(0.703);
39  add_FEMC_Cut(0.807);
40  add_FEMC_Cut(0.95);
41 
42  /* ---------------------------------------- */
43 
44  // Momentums //
45  vector< float > v_momenta;
46  v_momenta.push_back(1.0);
47  v_momenta.push_back(2.0);
48  v_momenta.push_back(5.0);
49  v_momenta.push_back(10.0);
50  v_momenta.push_back(20.0);
51 
52  // Marker Styles //
53  vector< int > v_styles;
54  v_styles.push_back(20);
55  v_styles.push_back(24);
56  v_styles.push_back(21);
57  v_styles.push_back(25);
58  v_styles.push_back(22);
59  v_styles.push_back(26);
60  v_styles.push_back(23);
61  v_styles.push_back(27);
62  v_styles.push_back(34);
63  v_styles.push_back(28);
64  v_styles.push_back(29);
65  v_styles.push_back(30);
66 
67  /* Loading Macros */
68 
69  gROOT->LoadMacro("efficiency_purity_analysisMacro.C");
70  gROOT->LoadMacro("/sphenix/user/gregtom3/SBU/research/macros/macros/sPHENIXStyle/sPhenixStyle.C");
72 
73  /* Useful local variables */
74  float v_rate_tp; //true positive rate value holder
75  float v_rate_fp; //false positive rate value holder
76 
77  //------------------------------------------------------------//
78  //-------- CEMC Plots --------//
79  //------------------------------------------------------------//
80 
81  //CEMC Graph Initialization//
82  const int size_C=2*v_cuts_CEMC.size(); //T.P & F.P graph per cut
83  TCanvas *c_CEMC = new TCanvas("cc","cc",1000,600);
84  TGraph *graphs_CEMC[size_C]; //Array of TGraph pointers
85  auto legend_CEMC_t = new TLegend(0.51,0.2,1.0,0.9); //Legend
86  auto legend_CEMC_f = new TLegend(0.51,0.2,1.0,0.9); //Legend
87 
88  legend_CEMC_t->SetTextSize(0.06);
89  legend_CEMC_f->SetTextSize(0.06);
90  // Create the CEMC TGraphs ( 2 per cut )//
91  for(unsigned idx_cut=0;idx_cut < v_cuts_CEMC.size(); idx_cut++)
92  /* loop over cuts */
93  {
94  // Index for each graph in TGraph vector //
95  int idx_tp=idx_cut*2;
96  int idx_fp=idx_cut*2+1;
97 
98  // Create new TGraphs_CEMC //
99  graphs_CEMC[idx_tp]=new TGraph(v_momenta.size());
100  graphs_CEMC[idx_fp]=new TGraph(v_momenta.size());
101 
102  // Clear true positive and false positive //
103  v_rate_tp = v_rate_fp = 0.0;
104 
105  // Set Points in TGraph for True Positive and False Positive per momenta//
106  for ( unsigned idx_p = 0; idx_p < v_momenta.size(); idx_p++ )
107  {
108  v_rate_tp=get_true_positive("C",v_cuts_CEMC.at(idx_cut),v_momenta.at(idx_p));
109  v_rate_fp=get_false_positive("C",v_cuts_CEMC.at(idx_cut),v_momenta.at(idx_p));
110  graphs_CEMC[idx_tp]->SetPoint(idx_p,v_momenta[idx_p],v_rate_tp);
111  graphs_CEMC[idx_fp]->SetPoint(idx_p,v_momenta[idx_p],v_rate_fp);
112  }
113 
114  // Tweak Graph Style //
115  graphs_CEMC[idx_tp]->SetMarkerColor(kGreen+2);
116  graphs_CEMC[idx_fp]->SetMarkerColor(kRed+2);
117 
118  graphs_CEMC[idx_tp]->SetLineColor(graphs_CEMC[idx_tp]->GetMarkerColor());
119  graphs_CEMC[idx_fp]->SetLineColor(graphs_CEMC[idx_fp]->GetMarkerColor());
120 
121  graphs_CEMC[idx_tp]->SetMarkerStyle(v_styles.at(idx_cut));
122  graphs_CEMC[idx_fp]->SetMarkerStyle(v_styles.at(idx_cut));
123 
124  graphs_CEMC[idx_tp]->SetMarkerSize(1); // Possibly too large? //
125  graphs_CEMC[idx_fp]->SetMarkerSize(1);
126 
127  graphs_CEMC[idx_tp]->SetTitle("");
128  graphs_CEMC[idx_fp]->SetTitle("");
129 
130 
131 
132  legend_CEMC_t->AddEntry(graphs_CEMC[idx_tp],v_cuts_CEMC_short.at(idx_cut)+" True Positive","P");
133  legend_CEMC_f->AddEntry(graphs_CEMC[idx_fp],v_cuts_CEMC_short.at(idx_cut)+ " False Positive","P");
134  }
135 
136  // Draw all graphs on same canvas (CEMC) //
137  c_CEMC->Divide(1,2);
138  for(int i = 0; i<size_C; i++)
139  {
140  if(i==0)
141  {
142  c_CEMC->cd(1);
143  //gPad->SetBottomMargin(0.2);
144  gPad->SetRightMargin(0.5);
145  graphs_CEMC[0]->Draw("PA");
146  graphs_CEMC[0]->GetXaxis()->SetTitle("True Momentum (GeV)");
147  graphs_CEMC[0]->GetYaxis()->SetRangeUser(0.7,1.0);
148  graphs_CEMC[0]->GetYaxis()->SetTitle("Rate");
149  }
150  else if(i==1)
151  {
152  c_CEMC->cd(2);
153  //gPad->SetTopMargin(0.2);
154  gPad->SetRightMargin(0.5);
155  graphs_CEMC[1]->GetXaxis()->SetTitle("True Momentum (GeV)");
156  graphs_CEMC[1]->GetYaxis()->SetRangeUser(-0.01,0.11);
157  graphs_CEMC[1]->GetYaxis()->SetTitle("Rate");
158  graphs_CEMC[1]->Draw("PA");
159  }
160  else if (i%2==0)
161  {
162  c_CEMC->cd(1);
163  graphs_CEMC[i]->Draw("PSame");
164  }
165  else
166  {
167  c_CEMC->cd(2);
168  graphs_CEMC[i]->Draw("P");
169  }
170 
171  }
172  c_CEMC->cd(1);
173  legend_CEMC_t->Draw();
174  c_CEMC->cd(2);
175  legend_CEMC_f->Draw();
176 
177  //------------------------------------------------------------//
178  //-------- EEMC Plots --------//
179  //------------------------------------------------------------//
180 
181 
182  // Create EEMC Graphs (2 per cut) //
183  const int size_E=v_cuts_EEMC.size(); //T.P & F.P graph per cut
184  TCanvas *c_EEMC = new TCanvas("ce","ce",1000,600);
185  TGraph *graphs_EEMC[size_E]; //Array of TGraph pointers
186  auto legend_EEMC_t = new TLegend(0.51,0.2,1,0.9); //Legend
187  auto legend_EEMC_f = new TLegend(0.51,0.2,1,0.9); //Legend
188 
189  legend_EEMC_t->SetTextSize(0.06);
190  legend_EEMC_f->SetTextSize(0.06);
191  for(unsigned idx_cut=0;idx_cut < v_cuts_EEMC.size()/2; idx_cut++)
192  // loop over cuts //
193  {
194  // Index for each graph in TGraph vector //
195  int idx_tp=idx_cut*2;
196  int idx_fp=idx_cut*2+1;
197 
198  // Create new TGraphs //
199  graphs_EEMC[idx_tp]=new TGraph(v_momenta.size());
200  graphs_EEMC[idx_fp]=new TGraph(v_momenta.size());
201 
202  // Clear true positive and false positive //
203  v_rate_tp = v_rate_fp = 0.0;
204 
205  // Fill data points for True Positve and False Positive Plot //
206  for ( unsigned idx_p = 0; idx_p < v_momenta.size(); idx_p++ )
207  {
208  TString cut = v_cuts_EEMC.at(idx_cut*2)+ getCut(v_momenta.at(idx_p)) + v_cuts_EEMC.at(idx_cut*2+1);
209  v_rate_tp=get_true_positive("E",cut,v_momenta.at(idx_p));
210  v_rate_fp=get_false_positive("E",cut,v_momenta.at(idx_p));
211  graphs_EEMC[idx_tp]->SetPoint(idx_p,v_momenta[idx_p],v_rate_tp);
212  graphs_EEMC[idx_fp]->SetPoint(idx_p,v_momenta[idx_p],v_rate_fp);
213  }
214 
215  graphs_EEMC[idx_tp]->SetMarkerColor(kGreen+2);
216  graphs_EEMC[idx_fp]->SetMarkerColor(kRed+2);
217 
218  graphs_EEMC[idx_tp]->SetLineColor(graphs_EEMC[idx_tp]->GetMarkerColor());
219  graphs_EEMC[idx_fp]->SetLineColor(graphs_EEMC[idx_fp]->GetMarkerColor());
220 
221  graphs_EEMC[idx_tp]->SetMarkerStyle(v_styles.at(idx_cut));
222  graphs_EEMC[idx_fp]->SetMarkerStyle(v_styles.at(idx_cut));
223 
224  graphs_EEMC[idx_tp]->SetMarkerSize(1); // Possibly too large? //
225  graphs_EEMC[idx_fp]->SetMarkerSize(1);
226 
227  graphs_EEMC[idx_tp]->SetTitle("");
228  graphs_EEMC[idx_fp]->SetTitle("");
229 
230  legend_EEMC_t->AddEntry(graphs_EEMC[idx_tp],v_cuts_EEMC_short.at(idx_cut)+" True Positive","P");
231  legend_EEMC_f->AddEntry(graphs_EEMC[idx_tp],v_cuts_EEMC_short.at(idx_cut)+" False Positive","P");
232  }
233  // Draw all graphs on same canvas (EEMC) //
234  c_EEMC->Divide(1,2);
235  for(int i = 0; i<size_E; i++)
236  {
237  if(i==0)
238  {
239  c_EEMC->cd(1);
240  gPad->SetRightMargin(0.5);
241  graphs_EEMC[0]->Draw("PA");
242  graphs_EEMC[0]->GetXaxis()->SetTitle("True Momentum (GeV)");
243  graphs_EEMC[0]->GetYaxis()->SetRangeUser(0.7,1.0);
244  graphs_EEMC[0]->GetYaxis()->SetTitle("Rate");
245  }
246  else if (i==1)
247  {
248  c_EEMC->cd(2);
249  gPad->SetRightMargin(0.5);
250  graphs_EEMC[1]->Draw("PA");
251  graphs_EEMC[1]->GetXaxis()->SetTitle("True Momentum (GeV)");
252  graphs_EEMC[1]->GetYaxis()->SetRangeUser(-0.01,0.11);
253  graphs_EEMC[1]->GetYaxis()->SetTitle("Rate");
254  }
255  else if (i%2==0)
256  {
257  c_EEMC->cd(1);
258  graphs_EEMC[i]->Draw("PSame");
259  }
260  else
261  {
262  c_EEMC->cd(2);
263  graphs_EEMC[i]->Draw("P");
264  }
265  }
266  c_EEMC->cd(1);
267  legend_EEMC_t->Draw();
268  c_EEMC->cd(2);
269  legend_EEMC_f->Draw();
270 
271  //------------------------------------------------------------//
272  //-------- FEMC Plots --------//
273  //------------------------------------------------------------//
274 
275 
276  // Create FEMC Graphs (2 per cut) //
277  const int size_F=v_cuts_FEMC.size(); //T.P & F.P graph per cut
278  TCanvas *c_FEMC = new TCanvas("cf","cf",1000,600);
279  TGraph *graphs_FEMC[size_F]; //Array of TGraph pointers
280  auto legend_FEMC_t = new TLegend(0.51,0.2,1.0,0.9); //Legend
281  auto legend_FEMC_f = new TLegend(0.51,0.2,1.0,0.9); //Legend
282 
283  legend_FEMC_t->SetTextSize(0.06);
284  legend_FEMC_f->SetTextSize(0.06);
285 
286  for(unsigned idx_cut=0;idx_cut < v_cuts_FEMC.size()/2; idx_cut++)
287  // loop over cuts //
288  {
289  // Index for each graph in TGraph vector //
290  int idx_tp=idx_cut*2;
291  int idx_fp=idx_cut*2+1;
292 
293  // Create new TGraphs //
294  graphs_FEMC[idx_tp]=new TGraph(v_momenta.size());
295  graphs_FEMC[idx_fp]=new TGraph(v_momenta.size());
296 
297  // Clear true positive and false positive //
298  v_rate_tp = v_rate_fp = 0.0;
299 
300  // Fill data points for True Positve and False Positive Plot //
301  for ( unsigned idx_p = 0; idx_p < v_momenta.size(); idx_p++ )
302  {
303  TString cut = v_cuts_FEMC.at(idx_cut*2)+ getCut(v_momenta.at(idx_p)) + v_cuts_FEMC.at(idx_cut*2+1);
304  v_rate_tp=get_true_positive("F",cut,v_momenta.at(idx_p));
305  v_rate_fp=get_false_positive("F",cut,v_momenta.at(idx_p));
306  graphs_FEMC[idx_tp]->SetPoint(idx_p,v_momenta[idx_p],v_rate_tp);
307  graphs_FEMC[idx_fp]->SetPoint(idx_p,v_momenta[idx_p],v_rate_fp);
308  }
309 
310  graphs_FEMC[idx_tp]->SetMarkerColor(kGreen+2);
311  graphs_FEMC[idx_fp]->SetMarkerColor(kRed+2);
312 
313  graphs_FEMC[idx_tp]->SetLineColor(graphs_FEMC[idx_tp]->GetMarkerColor());
314  graphs_FEMC[idx_fp]->SetLineColor(graphs_FEMC[idx_fp]->GetMarkerColor());
315 
316  graphs_FEMC[idx_tp]->SetMarkerStyle(v_styles.at(idx_cut));
317  graphs_FEMC[idx_fp]->SetMarkerStyle(v_styles.at(idx_cut));
318 
319  graphs_FEMC[idx_tp]->SetMarkerSize(1); // Possibly too large? //
320  graphs_FEMC[idx_fp]->SetMarkerSize(1);
321 
322  graphs_FEMC[idx_tp]->SetTitle("");
323  graphs_FEMC[idx_fp]->SetTitle("");
324 
325  legend_FEMC_t->AddEntry(graphs_FEMC[idx_tp],v_cuts_FEMC_short.at(idx_cut)+" True Positive","P");
326  legend_FEMC_f->AddEntry(graphs_FEMC[idx_fp],v_cuts_FEMC_short.at(idx_cut)+ " False Positive","P");
327  }
328 
329  // Draw all graphs on same canvas (FEMC) //
330  c_FEMC->Divide(1,2);
331  for(int i = 0; i<size_F; i++)
332  {
333  if(i==0)
334  {
335  c_FEMC->cd(1);
336  //gPad->SetBottomMargin(0.2);
337  gPad->SetRightMargin(0.5);
338  graphs_FEMC[0]->Draw("PA");
339  graphs_FEMC[0]->GetXaxis()->SetTitle("True Momentum (GeV)");
340  graphs_FEMC[0]->GetYaxis()->SetRangeUser(0.7,1.0);
341  graphs_FEMC[0]->GetYaxis()->SetTitle("Rate");
342  }
343  else if(i==1)
344  {
345  c_FEMC->cd(2);
346  //gPad->SetTopMargin(0.2);
347  gPad->SetRightMargin(0.5);
348  graphs_FEMC[1]->GetXaxis()->SetTitle("True Momentum (GeV)");
349  graphs_FEMC[1]->GetYaxis()->SetRangeUser(-0.01,0.11);
350  graphs_FEMC[1]->GetYaxis()->SetTitle("Rate");
351  graphs_FEMC[1]->Draw("PA");
352  }
353  else if (i%2==0)
354  {
355  c_FEMC->cd(1);
356  graphs_FEMC[i]->Draw("PSame");
357  }
358  else
359  {
360  c_FEMC->cd(2);
361  graphs_FEMC[i]->Draw("PSame");
362  }
363  }
364  c_FEMC->cd(1);
365  legend_FEMC_t->Draw();
366  c_FEMC->cd(2);
367  legend_FEMC_f->Draw();
368 
369  return 0;
370 }
371 
372 void add_CEMC_Cut(double ep, double prob)
373 {
374  //TStrings necessary to represent cuts//
375 
376  TString base_ep = "em_cluster_e/em_track_ptotal > ";
377  TString base_prob = "em_cluster_prob >= ";
378 
379  //TString necessary to print legible legends//
380 
381  TString base_ep_short = "E/P > ";
382  TString base_prob_short = "Prob >= ";
383  char str_1[40];
384  char str_2[40];
385  sprintf(str_1,"%.2f", ep);
386  sprintf(str_2,"%.3f", prob);
387 
388  //Test which cuts to use//
389 
390  bool do_ep=true; bool do_prob=true;
391  if(ep==0)
392  do_ep=false;
393  if(prob==0)
394  do_prob=false;
395 
396  //Add cuts to vector//
397 
398  if(do_prob&&!do_ep)
399  {
400  v_cuts_CEMC.push_back(TString(base_prob) += prob);
401  v_cuts_CEMC_short.push_back((TString(base_prob_short)+=str_2)+" True Positive");
402  }
403  else if(do_ep&&!do_prob)
404  {
405  v_cuts_CEMC.push_back(TString(base_ep) += ep);
406  v_cuts_CEMC_short.push_back((TString(base_ep_short)+=str_1)+" True Positive");
407  }
408  else
409  {
410  v_cuts_CEMC.push_back((TString(base_ep) += ep) + "&&" + (TString(base_prob) += prob));
411  v_cuts_CEMC_short.push_back((TString(base_ep_short) += str_1) + "&&" +(TString(base_prob_short) += str_2));
412  }
413 }
414 
415 void add_EEMC_Cut(double ep)
416 {
417  TString base_ep="em_cluster_e/";
418  TString base_ep_short = "E/P > ";
419  char str_1[40];
420  sprintf(str_1,"%.2f", ep);
421 
422  v_cuts_EEMC.push_back(base_ep); // 1st element : 'em_cluster_e/'
423  v_cuts_EEMC.push_back(TString(">")+= str_1); //: 'em_cluster_e/*GeV*>cut'
424 
425  v_cuts_EEMC_short.push_back(((TString(base_ep_short)) += str_1));
426 }
427 
428 void add_FEMC_Cut(double ep)
429 {
430  TString base_ep="em_cluster_e/";
431  TString base_ep_short = "E/P > ";
432  char str_1[40];
433  sprintf(str_1,"%.2f", ep);
434 
435  v_cuts_FEMC.push_back(base_ep); // 1st element : 'em_cluster_e/'
436  v_cuts_FEMC.push_back(TString(">")+= str_1); //: 'em_cluster_e/*GeV*>cut'
437 
438  v_cuts_FEMC_short.push_back(((TString(base_ep_short)) += str_1));
439 }