Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
makePlot_truth_LoI_Momentum.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file makePlot_truth_LoI_Momentum.C
1 
8 int makePlot_truth_LoI_Momentum(
9  const TString infile = "../../data/Sample_DISReco_ep.root",
10  //const TString infile = "../../data/EventGenAna_Pythia6_DIS_10x250_100k.root",
11  const bool save_figures = true
12  )
13 {
14  gStyle->SetOptStat(kFALSE);
15 
16  /*--------------Get Input File--------------*/
17 
18  TFile *f_truth = new TFile( infile, "OPEN");
19  TTree* t_truth = (TTree*)f_truth->Get("event_truth");
20  t_truth->Print();
21 
22  /*------------------------------------------*/
23  /*--------------Define Cuts-----------------*/
24 
25  //TCut Q2_cut = "evtgen_Q2 > 1";
26 
27  TCut electron_cut = "em_evtgen_pid == 11";
28  TCut hadron_cut = "abs(em_evtgen_pid) > 100";
29  TCut proton_cut = "em_evtgen_pid == 2212";
30  TCut neutron_cut = "em_evtgen_pid == 2112";
31  TCut Kaon_cut = "abs(em_evtgen_pid) == 321 || em_evtgen_pid == 311";
32  TCut Pion_charged_cut = "abs(em_evtgen_pid) == 211";
33  //TCut Pion_cut = "abs(em_evtgen_pid) == 211 || em_evtgen_pid == 111";
34  TCut photon_cut = "em_evtgen_pid == 22";
35 
36  TCut mother_cut = "1"; //"p.fParent == 0";
37  TCut status_cut = "1"; //"p.fKS < 10";
38 
39  TCut eta_cut_n3n2 = "(-1 * log( tan( em_evtgen_theta / 2.0 ) )) > -3 && (-1 * log( tan( em_evtgen_theta / 2.0 ) )) < -2";
40  TCut eta_cut_n2n1 = "(-1 * log( tan( em_evtgen_theta / 2.0 ) )) > -2 && (-1 * log( tan( em_evtgen_theta / 2.0 ) )) < -1";
41  TCut eta_cut_n1z0 = "(-1 * log( tan( em_evtgen_theta / 2.0 ) )) > -1 && (-1 * log( tan( em_evtgen_theta / 2.0 ) )) < -0";
42 
43  /*------------------------------------------*/
44  /*-------Momentum vs. Pseudorapidity--------*/
45  /*------------------------------------------*/
46 
47  /*------------------------------------------*/
48  /*---------Electrons (LoI Fig. 2-1)---------*/
49  /*------------------------------------------*/
50  TH2F *h_peta_e = new TH2F("h_peta_e", "", 100, -4, 3, 100, 0, 50 ); //250x015 10M
51 
52  TCanvas *c_peta_e = new TCanvas( "c_peta_e" );
53  c_peta_e->SetLogz();
54 
55  t_truth->Draw("em_evtgen_ptotal:(-1 * log( tan( em_evtgen_theta / 2.0 ) ))>>h_peta_e", electron_cut && mother_cut && "evtgen_Q2 > 1", "colz");
56 
57  h_peta_e->GetXaxis()->SetTitle("Pseudorapidity #eta");
58  h_peta_e->GetYaxis()->SetTitle("Electron Momentum p_{e-} [GeV]");
59 
60  if ( save_figures )
61  {
62  c_peta_e->Print("plots/Pythia_peta_e.eps");
63  c_peta_e->Print("plots/Pythia_peta_e.png");
64  }
65 
66  /*------------------------------------------*/
67  /*----------Hadrons (LoI Fig. 2-4)----------*/
68  /*------------------------------------------*/
69  TH2F *h_p_eta_h = new TH2F("h_p_eta_h", "", 60,-6,6, 200,0,225); //250x015 10M
70 
71  TCanvas *c_p_eta_h = new TCanvas( "c_p_eta_h" );
72  t_truth->Draw("em_evtgen_ptotal:(-1 * log( tan( em_evtgen_theta / 2.0 ) ))>>h_p_eta_h",hadron_cut && "evtgen_Q2 > 1 && 0.01 < evtgen_y < 0.80", "colz");
73 
74  h_p_eta_h->GetXaxis()->SetTitle("Pseudorapidity #eta");
75  h_p_eta_h->GetYaxis()->SetTitle("Hadron Momentum p_{Hadron} [GeV]");
76  c_p_eta_h->SetLogz();
77 
78  if ( save_figures )
79  {
80  c_p_eta_h->Print("plots/Pythia_peta_h.eps");
81  c_p_eta_h->Print("plots/Pythia_peta_h.png");
82  }
83 
84  /*------------Momentum Spectra (Figure 2.2)--------------*/
85 
86  /*-------- -3 < Eta < -2 ---------*/
87 
88  TH1F* hp_e_n3n2 = new TH1F("hp_e_n3n2", "dN/dp vs. p", 60, 0, 30);
89  TH1F* hp_p_n3n2 = new TH1F("hp_p_n3n2", "dN/dp vs. p", 60, 0, 30);
90  TH1F* hp_y_n3n2 = new TH1F("hp_y_n3n2", "dN/dp vs. p", 60, 0, 30);
91 
92  TCanvas *cp_e_n3n2 = new TCanvas("cp_e_n3n2");
93  t_truth->Draw("em_evtgen_ptotal>>hp_e_n3n2", electron_cut && eta_cut_n3n2 && "evtgen_Q2 > 0.01", "goff");
94  t_truth->Draw("em_evtgen_ptotal>>hp_p_n3n2", Pion_charged_cut && eta_cut_n3n2 && "evtgen_Q2 > 0.01", "goff");
95  t_truth->Draw("em_evtgen_ptotal>>hp_y_n3n2", photon_cut && eta_cut_n3n2 && status_cut && "evtgen_Q2 > 0.01", "goff");
96 
97 
98  TH1F* htmp_n3n2 = hp_e_n3n2->Clone();
99  htmp_n3n2->SetTitle("");
100  htmp_n3n2->GetXaxis()->SetTitle("p [GeV]");
101  htmp_n3n2->GetYaxis()->SetTitle("dN/dp");
102  htmp_n3n2->SetMaximum( 0.99e7);
103  htmp_n3n2->Draw();
104 
105  hp_e_n3n2->SetLineColor(2);
106  hp_e_n3n2->Draw("same");
107  hp_p_n3n2->SetLineColor(1);
108  hp_p_n3n2->Draw("same");
109  hp_y_n3n2->SetLineColor(4);
110  hp_y_n3n2->Draw("same");
111  cp_e_n3n2->SetLogy();
112 
113  TLegend* leg_n3n2 = new TLegend(0.53,0.67,0.73,0.90);
114  hp_e_n3n2->SetTitle("DIS electron");
115  hp_p_n3n2->SetTitle("#pi#pm");
116  hp_y_n3n2->SetTitle("Photons");
117  leg_n3n2->AddEntry(hp_e_n3n2, "", "L");
118  leg_n3n2->AddEntry(hp_p_n3n2, "", "L");
119  leg_n3n2->AddEntry(hp_y_n3n2, "", "L");
120  leg_n3n2->SetTextSize(0.04);
121  leg_n3n2->Draw();
122 
123  if ( save_figures )
124  {
125  cp_e_n3n2->Print("plots/Pythia_MomSpec_n3n2.eps");
126  cp_e_n3n2->Print("plots/Pythia_MomSpec_n3n2.png");
127  }
128 
129  /*Temporary Electron Purity Test*/
130 
131 
132  /*-------- -2 < Eta < -1 ---------*/
133 
134  TH1F* hp_e_n2n1 = new TH1F("hp_e_n2n1", "dN/dp vs. p", 60, 0, 30);
135  TH1F* hp_p_n2n1 = new TH1F("hp_p_n2n1", "dN/dp vs. p", 60, 0, 30);
136  TH1F* hp_y_n2n1 = new TH1F("hp_y_n2n1", "dN/dp vs. p", 60, 0, 30);
137 
138  TCanvas *cp_e_n2n1 = new TCanvas("cp_e_n2n1");
139  t_truth->Draw("em_evtgen_ptotal>>hp_e_n2n1", electron_cut && eta_cut_n2n1 && "evtgen_Q2 > 0.01", "goff");
140  t_truth->Draw("em_evtgen_ptotal>>hp_p_n2n1", Pion_charged_cut && eta_cut_n2n1 && "evtgen_Q2 > 0.01", "goff");
141  t_truth->Draw("em_evtgen_ptotal>>hp_y_n2n1", photon_cut && eta_cut_n2n1 && status_cut && "evtgen_Q2 > 0.01", "goff");
142 
143  TH1F* htmp_n2n1 = hp_e_n2n1->Clone();
144  htmp_n2n1->SetTitle("");
145  htmp_n2n1->GetXaxis()->SetTitle("p [GeV]");
146  htmp_n2n1->GetYaxis()->SetTitle("dN/dp");
147  htmp_n2n1->SetMaximum( 0.99e7);
148  htmp_n2n1->Draw();
149 
150  hp_e_n2n1->SetLineColor(2);
151  hp_e_n2n1->Draw("same");
152  hp_p_n2n1->SetLineColor(1);
153  hp_p_n2n1->Draw("same");
154  hp_y_n2n1->SetLineColor(4);
155  hp_y_n2n1->Draw("same");
156  cp_e_n2n1->SetLogy();
157 
158  TLegend* leg_n2n1 = new TLegend(0.53,0.67,0.73,0.90);
159  hp_e_n2n1->SetTitle("DIS electron");
160  hp_p_n2n1->SetTitle("#pi#pm");
161  hp_y_n2n1->SetTitle("Photons");
162  leg_n2n1->AddEntry(hp_e_n2n1, "", "L");
163  leg_n2n1->AddEntry(hp_p_n2n1, "", "L");
164  leg_n2n1->AddEntry(hp_y_n2n1, "", "L");
165  leg_n2n1->SetTextSize(0.04);
166  leg_n2n1->Draw();
167 
168  if ( save_figures )
169  {
170  cp_e_n2n1->Print("plots/Pythia_MomSpec_n2n1.eps");
171  cp_e_n2n1->Print("plots/Pythia_MomSpec_n2n1.png");
172  }
173 
174  /*-------- -1 < Eta < 0 ---------*/
175 
176  TH1F* hp_e_n1z0 = new TH1F("hp_e_n1z0", "dN/dp vs. p", 60, 0, 30);
177  TH1F* hp_p_n1z0 = new TH1F("hp_p_n1z0", "dN/dp vs. p", 60, 0, 30);
178  TH1F* hp_y_n1z0 = new TH1F("hp_y_n1z0", "dN/dp vs. p", 60, 0, 30);
179 
180  TCanvas *cp_e_n1z0 = new TCanvas("cp_e_n1z0");
181  t_truth->Draw("em_evtgen_ptotal>>hp_e_n1z0", electron_cut && eta_cut_n1z0 && "evtgen_Q2 > 0.01", "goff");
182  t_truth->Draw("em_evtgen_ptotal>>hp_p_n1z0", Pion_charged_cut && eta_cut_n1z0 && "evtgen_Q2 > 0.01", "goff");
183  t_truth->Draw("em_evtgen_ptotal>>hp_y_n1z0", photon_cut && eta_cut_n1z0 && status_cut && "evtgen_Q2 > 0.01", "goff");
184 
185  TH1F* htmp_n1z0 = hp_e_n1z0->Clone();
186  htmp_n1z0->SetTitle("");
187  htmp_n1z0->GetXaxis()->SetTitle("p [GeV]");
188  htmp_n1z0->GetYaxis()->SetTitle("dN/dp");
189  htmp_n1z0->SetMaximum( 0.99e7 );
190  htmp_n1z0->Draw();
191 
192  hp_e_n1z0->SetLineColor(2);
193  hp_e_n1z0->Draw("same");
194  hp_p_n1z0->SetLineColor(1);
195  hp_p_n1z0->Draw("same");
196  hp_y_n1z0->SetLineColor(4);
197  hp_y_n1z0->Draw("same");
198  cp_e_n1z0->SetLogy();
199 
200  TLegend* leg_n1z0 = new TLegend(0.53,0.67,0.73,0.90);
201  hp_e_n1z0->SetTitle("DIS electron");
202  hp_p_n1z0->SetTitle("#pi#pm");
203  hp_y_n1z0->SetTitle("Photons");
204  leg_n1z0->AddEntry(hp_e_n1z0, "", "L");
205  leg_n1z0->AddEntry(hp_p_n1z0, "", "L");
206  leg_n1z0->AddEntry(hp_y_n1z0, "", "L");
207  leg_n1z0->SetTextSize(0.04);
208  leg_n1z0->Draw();
209 
210  if ( save_figures )
211  {
212  cp_e_n1z0->Print("plots/Pythia_MomSpec_n1z0.eps");
213  cp_e_n1z0->Print("plots/Pythia_MomSpec_n1z0.png");
214  }
215 
216  return 0 ;
217 }