source: trunk/FACT++/spectrum/display.C@ 20041

Last change on this file since 20041 was 19971, checked in by tbretz, 5 years ago
Added a binning in Impact Parameter for crosschecks
File size: 9.7 KB
Line 
1void display(const char *filename="spectrum.hist.root")
2{
3 gEnv->SetValue("Canvas.ShowEventStatus", 1);
4 gEnv->SetValue("TFile.Recover", 0);
5
6 // --------------------------------
7
8 TFile file(filename);
9 if (file.IsZombie())
10 return;
11
12 // --------------------------------
13
14 cout << endl;
15
16 file.ReadAll();
17 file.ls();
18
19 cout << endl;
20
21 // --------------------------------
22
23 TH1D *h1 = 0;
24 TH2D *h2 = 0;
25
26 TCanvas *c = 0;
27
28 gStyle->SetErrorX(0.5);
29
30 TF1 fx("fx", "x", 0, 90);
31 fx.SetLineStyle(kDashed);
32 fx.SetLineColor(kBlack);
33
34 TF1 f0("f0", "[0]", 0, 90);
35 f0.SetParameter(0, 0);
36 f0.SetLineStyle(kDashed);
37 f0.SetLineColor(kBlack);
38
39 TF1 f("spectrum", "[0]*pow(pow(10, x)/1000, -[1])", 2.8, 4.5);
40 f.SetLineColor(kBlue);
41
42 // --------------------------------
43
44 c = new TCanvas("Zd Weights", "Zd Weights");
45 c->Divide(2,2);
46 c->cd(1);
47 file.GetObject("Zd/OnTime", h1);
48 h1->SetTitle("OnTime (hist), Scaled excess (markers)");
49 double time = h1->Integral();
50 h1->DrawCopy();
51 file.GetObject("Data/Theta/Excess", h1);
52 h1->Scale(time/h1->Integral());
53 h1->DrawCopy("same");
54 c->cd(2);
55 file.GetObject("Zd/CountN", h1);
56 h1->DrawCopy("P");
57 c->cd(3);
58 file.GetObject("Zd/ZdWeight", h1);
59 h1->DrawCopy("P");
60
61 // --------------------------------
62
63 c = new TCanvas("Zenith Angle", "Zenith Angle");
64 c->Divide(2,2);
65
66 c->cd(1);
67 gPad->SetGridy();
68 file.GetObject("MC/theta/BiasW", h1);
69 h1->DrawCopy("P");
70 f0.DrawCopy("same");
71 c->cd(2);
72 gPad->SetGridy();
73 file.GetObject("MC/theta/ResolutionW", h1);
74 h1->DrawCopy("][ P");
75 f0.DrawCopy("same");
76
77 c->cd(3);
78 gPad->SetGridy();
79 file.GetObject("MC/theta/TriggerEfficiencyW", h1);
80 h1->SetTitle("Trigger (black) and Cut Efficiency (blue)");
81 h1->DrawCopy("P");
82 file.GetObject("MC/theta/CutEfficiencyW", h1);
83 h1->SetMarkerColor(kBlue);
84 h1->DrawCopy("P same");
85
86 c->cd(4);
87 gPad->SetLogy();
88 file.GetObject("Data/Theta/Spectrum", h1);
89 h1->DrawCopy("P");
90 file.GetObject("Data/Theta/RolkeUL", h1);
91 h1->SetMarkerStyle(23);
92 h1->DrawCopy("P same");
93
94 // --------------------------------
95
96 c = new TCanvas("Impact", "Impact");
97 c->Divide(2,2);
98 c->cd(1);
99 file.GetObject("MC/dense/TrueEnergy/Impact", h2);
100 h2->DrawCopy("colz");
101 h2->DrawCopy("same");
102 c->cd(2);
103 file.GetObject("MC/sparse/TrueEnergy/Impact", h2);
104 h2->DrawCopy("colz");
105 h2->DrawCopy("same");
106 c->cd(3);
107 file.GetObject("MC/theta/Impact", h2);
108 h2->DrawCopy("colz");
109 h2->DrawCopy("same");
110
111 // --------------------------------
112
113 c = new TCanvas("Migration", "Migration");
114 c->Divide(2,1);
115 c->cd(1);
116 file.GetObject("MC/dense/Migration", h2);
117 h2->DrawCopy("colz");
118 fx.DrawCopy("same");
119 c->cd(2);
120 file.GetObject("MC/sparse/Migration", h2);
121 h2->DrawCopy("colz");
122 fx.DrawCopy("same");
123
124 // --------------------------------
125 c = new TCanvas("Resolution", "Resolution");
126 c->Divide(2,2);
127 c->cd(1);
128 gPad->SetGridy();
129 file.GetObject("MC/dense/TrueEnergy/BiasW", h1);
130 h1->DrawCopy("P");
131 f0.DrawCopy("same");
132 c->cd(2);
133 gPad->SetGridy();
134 file.GetObject("MC/dense/EstimatedEnergy/BiasW", h1);
135 h1->DrawCopy("P");
136 f0.DrawCopy("same");
137 c->cd(3);
138 file.GetObject("MC/dense/TrueEnergy/ResolutionW", h1);
139 h1->DrawCopy("][ P");
140 c->cd(4);
141 file.GetObject("MC/dense/EstimatedEnergy/ResolutionW", h1);
142 h1->DrawCopy("][ P");
143
144 // --------------------------------
145
146 c = new TCanvas("Efficiency", "Efficiency");
147 c->Divide(2,2);
148 c->cd(1);
149 file.GetObject("MC/dense/TrueEnergy/TriggerEfficiencyW", h1);
150 h1->SetTitle("Trigger (black) and Cut (blue) efficiency");
151 h1->DrawCopy("P X0");
152 h1->DrawCopy("C hist same");
153 file.GetObject("MC/dense/TrueEnergy/CutEfficiencyW", h1);
154 h1->SetLineColor(kBlue);
155 h1->DrawCopy("P X0 same");
156 h1->DrawCopy("C hist same");
157
158 c->cd(2);
159 gPad->SetLogy();
160 file.GetObject("MC/dense/TrueEnergy/EffectiveAreaW", h1);
161 h1->SetTitle("Effective Collection Area");
162 h1->DrawCopy("P X0");
163 h1->DrawCopy("C hist same");
164
165 c->cd(3);
166
167 gPad->SetLogy();
168 file.GetObject("MC/dense/TrueEnergy/SimFluxW", h1);
169 h1->SetTitle("Simulated (blk), Triggered (blk), Cut (blu), Bg (blu), Reconstructed (red) spectrum");
170 h1->DrawCopy("P X0");
171 h1->DrawCopy("C hist same");
172
173 gPad->SetLogy();
174 file.GetObject("MC/dense/TrueEnergy/TrigFluxW", h1);
175 h1->DrawCopy("P X0 same");
176 h1->DrawCopy("C hist same");
177
178 gPad->SetLogy();
179 file.GetObject("MC/dense/TrueEnergy/ExcessFluxW", h1);
180 h1->DrawCopy("P X0 same");
181 h1->DrawCopy("C hist same");
182
183 TSpline3 spline(h1);
184 spline.SetLineColor(kBlue);
185 spline.DrawClone("same");
186
187 float max = 0;
188 float xmax = 0;
189 for (float x=2.5; x<5; x+=0.01)
190 {
191 if (spline.Eval(x)>max)
192 {
193 xmax = x;
194 max = spline.Eval(x);
195 }
196 }
197
198 TMarker m;
199 m.SetMarkerColor(kBlue);
200 m.SetMarkerStyle(kStar);
201 m.DrawMarker(xmax, max);
202
203 cout << endl;
204 cout << "Threshold Max: " << pow(10, xmax) << endl;
205 cout << endl;
206
207 file.GetObject("MC/dense/EstimatedEnergy/ExcessFluxW", h1);
208 h1->SetLineColor(kRed);
209 h1->DrawCopy("X0 P same");
210 h1->DrawCopy("C hist same");
211
212 file.GetObject("MC/dense/TrueEnergy/BackgroundFluxW", h1);
213 h1->SetLineColor(kBlue);
214 h1->DrawCopy("X0 P same");
215 h1->DrawCopy("C hist same");
216
217 c->cd(4);
218 gPad->SetLogy();
219 file.GetObject("Data/Energy/Differential/Excess", h1);
220 h1->SetTitle("Measured Signal (blue) / Simulated Excess (black)");
221 h1->SetLineColor(kBlue);
222 h1->DrawCopy("P X0");
223 h1->DrawCopy("C hist same");
224 double scale = h1->Integral();
225 file.GetObject("MC/sparse/EstimatedEnergy/ExcessFluxW", h1);
226 h1->Scale(scale/h1->Integral());
227 h1->DrawCopy("X0 P same");
228 h1->DrawCopy("C hist same");
229
230 // --------------------------------
231 f0.SetLineColor(kRed);
232 // --------------------------------
233
234 c = new TCanvas("Integral Spectrum", "Integral Spectrum");
235 c->Divide(2,2);
236 c->cd(4);
237 gPad->SetLogy();
238 file.GetObject("Data/Energy/Differential/IntegratedSpectrum", h1);
239 h1->SetLineColor(kGray);
240 h1->SetMarkerColor(kGray);
241 if (h1->GetMinimum()<=0)
242 {
243 h1->SetMinimum(h1->GetMinimum(0)*0.1);
244 cout << "WARNING: Integral Spectrum contains negative flux points!\n" << endl;
245 }
246 h1->DrawCopy("P");
247 file.GetObject("Data/Energy/Integral/Spectrum", h1);
248 h1->DrawCopy("P same");
249
250 f.SetParameters(h1->GetMaximum(), 1.4);
251 h1->Fit(&f, "N0EM", "");
252 f.DrawCopy("same");
253
254 cout << endl;
255 cout << "ChiSq " << f.GetChisquare() << " / " << f.GetNDF() << endl;
256 cout << "Prob. " << f.GetProb() << endl;
257 cout << endl;
258
259 file.GetObject("Data/Energy/Integral/RolkeUL", h1);
260 h1->SetMarkerStyle(23);
261 h1->DrawCopy("P same");
262
263 //file.GetObject("Data/Energy/Integral/RolkeLL", h1);
264 //h1->SetMarkerStyle(22);
265 //h1->DrawCopy("P same");
266
267 c->cd(1);
268 gPad->SetGridx();
269 file.GetObject("Data/Energy/Integral/SigmaFlux", h1);
270 h1->SetMaximum(5);
271 h1->SetMarkerStyle(kFullDotLarge);
272 h1->DrawCopy("P");
273 f0.SetParameter(0, 1);
274 f0.DrawCopy("same");
275
276 c->cd(3);
277 gPad->SetGridx();
278 file.GetObject("Data/Energy/Integral/BackgroundI", h1);
279 h1->SetTitle("Integral Signal (black), Background (red) and Excess (blue)");
280 h1->Scale(5);
281 h1->SetMaximum(30);
282 h1->SetMarkerColor(kRed);
283 h1->SetLineColor(kRed);
284 h1->DrawCopy("P");
285 file.GetObject("Data/Energy/Integral/SignalI", h1);
286 h1->DrawCopy("P same");
287 file.GetObject("Data/Energy/Integral/ExcessI", h1);
288 h1->SetMarkerColor(kBlue);
289 h1->SetLineColor(kBlue);
290 h1->DrawCopy("P same");
291
292 f0.SetParameter(0, 10);
293 f0.DrawCopy("same");
294
295 // --------------------------------
296
297 c = new TCanvas("Differential Spectrum", "Differential Spectrum");
298 c->Divide(2,2);
299 c->cd(4);
300 gPad->SetLogy();
301 file.GetObject("Data/Energy/Differential/Spectrum", h1);
302 if (h1->GetMinimum()<=0)
303 {
304 h1->SetMinimum(h1->GetMinimum(0)*0.1);
305 cout << "WARNING: Differential Spectrum contains negative flux points!\n" << endl;
306 }
307 h1->DrawCopy("P");
308
309 f.SetParameters(h1->GetMaximum(), 2.4);
310 h1->Fit(&f, "N0EM", "");
311 f.DrawCopy("same");
312
313 cout << endl;
314 cout << "ChiSq " << f.GetChisquare() << " / " << f.GetNDF() << endl;
315 cout << "Prob. " << f.GetProb() << endl;
316 cout << endl;
317
318 file.GetObject("Data/Energy/Differential/RolkeUL", h1);
319 h1->SetMarkerStyle(23);
320 h1->DrawCopy("P same");
321
322 //file.GetObject("Data/Energy/Differential/RolkeLL", h1);
323 //h1->SetMarkerStyle(22);
324 //h1->DrawCopy("P same");
325
326 //file.GetObject("Data/Energy/Differential/FeldmanCousins", h1);
327 //h1->SetMarkerColor(kBlue);
328 //h1->SetMarkerStyle(23);
329 //h1->DrawCopy("P same");
330
331 c->cd(1);
332 gPad->SetGridx();
333 file.GetObject("Data/Energy/Differential/SigmaFlux", h1);
334 h1->SetMaximum(5);
335 h1->SetMarkerStyle(kFullDotLarge);
336 h1->DrawCopy("P");
337 f0.SetParameter(0, 1);
338 f0.DrawCopy("same");
339
340 c->cd(3);
341 gPad->SetGridx();
342 file.GetObject("Data/Energy/Differential/Background", h1);
343 h1->SetTitle("Differential Signal (black), Background (red) and Excess (blue)");
344 h1->Scale(5);
345 h1->SetMaximum(30);
346 h1->SetMarkerColor(kRed);
347 h1->SetLineColor(kRed);
348 h1->DrawCopy("P");
349 file.GetObject("Data/Energy/Differential/Signal", h1);
350 h1->DrawCopy("P same");
351 file.GetObject("Data/Energy/Differential/Excess", h1);
352 h1->SetMarkerColor(kBlue);
353 h1->SetLineColor(kBlue);
354 h1->DrawCopy("P same");
355
356 f0.SetParameter(0, 10);
357 f0.DrawCopy("same");
358}
Note: See TracBrowser for help on using the repository browser.