FairRoot/PandaRoot
crosstag.C
Go to the documentation of this file.
1 #include <algorithm>
2 #include <map>
3 #include <utility>
4 #include <iostream>
5 
6 #include "TFile.h"
7 #include "TTree.h"
8 #include "TString.h"
9 #include "TCanvas.h"
10 #include "TH1F.h"
11 #include "TH2F.h"
12 #include "TROOT.h"
13 #include "TEventList.h"
14 #include "TDirectory.h"
15 #include "TRegexp.h"
16 #include "TStyle.h"
17 #include "TColor.h"
18 #include "TPaletteAxis.h"
19 #include "TLegend.h"
20 
21 std::map<int,int> codeidx, drawline;
22 
23 Int_t tagline[60], mode, recmode, tag, tagall, tagm, codes[60], nlines;
24 
25 // -----------------------------------------------------------------------
26 
27 void palette2()
28 {
29  Double_t blu[3] = {1,1,1};
30  Double_t red[3] = {1,0.5,0};
31  Double_t grn[3] = {1,0.5,0};
32  Double_t stp[3] = {0,0.5,1};
33 
34  TColor::CreateGradientColorTable(3,stp,red,grn,blu,100);
35 }
36 
37 // -----------------------------------------------------------------------
38 
39 void setStyle()
40 {
41  gStyle->SetPadTopMargin(0.05);
42  gStyle->SetPadBottomMargin(0.1);
43  gStyle->SetPadLeftMargin(0.1);
44  gStyle->SetPadRightMargin(0.07);
45  gStyle->SetTitleH(0.035);
46  //gStyle->SetTitleX(0.1);
47  //gStyle->SetTitleW(0.9);
48  gStyle->SetPalette(1);
49 }
50 
51 // -----------------------------------------------------------------------
52 
53 void resizePalette(TH1* h)
54 {
55  TPaletteAxis *palette = (TPaletteAxis*)h->GetListOfFunctions()->FindObject("palette");
56  palette->SetX1NDC(0.94);
57  palette->SetX2NDC(0.96);
58  palette->SetY1NDC(0.1);
59  palette->SetY2NDC(0.95);
60 }
61 
62 // -----------------------------------------------------------------------
63 
64 void init(TTree *t, double sqs)
65 {
66  // prepare tree by setting branches
67  t->SetBranchAddress("mode",&mode);
68  t->SetBranchAddress("recmode",&recmode);
69  t->SetBranchAddress("tag",&tag);
70  t->SetBranchAddress("tagall",&tagall);
71  t->SetBranchAddress("tagm",&tagm);
72 
73  t->SetBranchStatus("*",0);
74  t->SetBranchStatus("tag*",1);
75  t->SetBranchStatus("mode",1);
76  t->SetBranchStatus("recmode",1);
77  t->SetBranchStatus("tag",1);
78  t->SetBranchStatus("tagall",1);
79  t->SetBranchStatus("tagm",1);
80 
81  TObjArray* branches = t->GetListOfBranches();
82 
83  nlines = 0;
84  TRegexp reg("tag[0-9][0-9][0-9]");
85 
86  for(int i=0; i<=branches->GetLast(); ++i)
87  {
88  TBranch* branch = (TBranch*)branches->UncheckedAt(i);
89  TString v=branch->GetName();
90  if (v(reg)!="")
91  {
92  int m = TString(v(3,3)).Atoi();
93  if (sqs!=3.0 || m>=400)
94  {
95  t->SetBranchAddress(v,&(tagline[nlines]));
96  //t->SetBranchStatus(v,1);
97  codeidx[m]=nlines;
98  codes[nlines] = m;
99 
100  nlines++;
101  }
102  }
103  }
104  codeidx[900]=nlines;
105  codes[nlines]=900;
106  nlines++;
107 
108  drawline[110] = 1;
109  drawline[120] = 1;
110  drawline[130] = 1;
111  drawline[140] = 1;
112  drawline[200] = 1;
113  drawline[400] = 1;
114  drawline[500] = 1;
115  drawline[600] = 1;
116 }
117 
118 // -----------------------------------------------------------------------
119 
120 void config_pad(TVirtualPad* p, double b=0.1, double r=0.07, double t=0.05, double l=0.1)
121 {
122  p->SetTopMargin(t);
123  p->SetBottomMargin(b);
124  p->SetLeftMargin(l);
125  p->SetRightMargin(r);
126  p->SetGridx();
127  p->SetGridy();
128 }
129 
130 // -----------------------------------------------------------------------
131 
132 void config_histo(TH1* h, double labs=0.018, TString titley="", TString titlex="")
133 {
134  h->GetXaxis()->SetLabelSize(labs);
135  h->GetXaxis()->SetLabelOffset(0.01);
136  h->GetXaxis()->SetTitleSize(0.03);
137  h->GetXaxis()->SetTitleOffset(1.5);
138 
139  h->GetYaxis()->SetLabelSize(labs);
140  h->GetYaxis()->SetLabelOffset(0.01);
141  h->GetYaxis()->SetTitleSize(0.03);
142  h->GetYaxis()->SetTitleOffset(1.5);
143 
144  h->GetZaxis()->SetLabelSize(0.018);
145  h->GetZaxis()->SetLabelOffset(0.005);
146 
147  h->SetStats(0);
148  h->SetContour(100);
149 
150  if (titley!="") h->SetYTitle(titley);
151  if (titlex!="") h->SetXTitle(titlex);
152 
153  //h->Sumw2();
154  h->SetLineWidth(2);
155 }
156 
157 // -----------------------------------------------------------------------
158 
159 void config_histo2d(TH1* h, TString titley="", TString titlex="")
160 {
161  double labs = 0.018+0.008*(57.-nlines)/43.;
162 
163  config_histo(h, labs, titley, titlex);
164 
165  for (int i=0;i<nlines-1;++i)
166  {
167  h->GetXaxis()->SetBinLabel(i+1,TString::Format("T%d",codes[i]));
168  h->GetYaxis()->SetBinLabel(nlines-i,TString::Format("M%d",codes[i]));
169  }
170 
171  h->GetYaxis()->SetBinLabel(1,"DPM");
172  h->GetXaxis()->SetBinLabel(nlines,"any");
173  h->GetXaxis()->LabelsOption("v");
174 }
175 
176 // -----------------------------------------------------------------------
177 
178 void config_histo1d(TH1* h, TString titley="", TString titlex="", int lincol=1, int fillcol=0, int fillstyle=0, TString label="M")
179 {
180  double labs = 0.04+0.01*(57.-nlines)/43.;
181 
182  config_histo(h, labs, titley, titlex);
183 
184  h->GetYaxis()->SetTitleSize(0.05);
185  h->GetYaxis()->SetTitleOffset(0.7);
186  h->GetYaxis()->SetLabelSize(0.04);
187 
188  h->GetXaxis()->SetTitleSize(0.045);
189 
190  for (int i=0;i<nlines-1;++i)
191  h->GetXaxis()->SetBinLabel(i+1,TString::Format("%s%d",label.Data(),codes[i]));
192 
193  h->GetXaxis()->SetBinLabel(nlines,"DPM");
194  h->GetXaxis()->LabelsOption("v");
195 
196  h->SetLineColor(lincol);
197  h->SetFillColor(fillcol);
198  h->SetFillStyle(fillstyle);
199 }
200 
201 // -----------------------------------------------------------------------
202 
203 int crosstag(TString fname, int fact=1, int saveplots=0, double sqs=-1.)
204 {
205  setStyle();
206 
207  TFile *f = new TFile(fname,"READ");
208  TTree *t = (TTree*)f->Get("ntpev");
209 
210  if (sqs<0.) sqs = TString(fname(fname.Index("/M")+2,3)).Atoi()/100.; // assume file name Mxxx_... with xxx = sqs*100
211  if (sqs<0. || sqs>6.) return;
212 
213 
214  // ********
215  // init branches in trees, global vars etc
216  // ********
217  init(t, sqs);
218 
219  TString sqsstr = TString::Format("%03d",int(sqs*100.));
220 
221  TCanvas *c1=new TCanvas("c1","c1",10,10,1000,900);
222  TCanvas *c3=new TCanvas("c3","c3",100,100,1000,900);
223  TCanvas *c2=new TCanvas("c2","c2",400,20,1000,900);
224 
225  c2->Divide(1,2);
226 
227  TFile *ff=0;
228  if (saveplots) ff=new TFile("crosstag_histos.root","UPDATE");
229 
230 
231  // ********
232  // 2D Plots
233  // ********
234 
235  // normalization histo
236  TH2F *hall = new TH2F("hall"+sqsstr,"normalization",nlines,0,nlines,nlines,0,nlines);
237 
238  // matrix of tag line vs data mode
239  TH2F *htag = new TH2F("htag"+sqsstr,TString::Format("Cross tags @ %.1f GeV",sqs),nlines,0,nlines,nlines,0,nlines);
240 
241  // matrix of tag line vs data mode, containing the exclusive triggers for certain data modes
242  TH2F *htagex = new TH2F("htagex"+sqsstr,TString::Format("Excl. cross tags @ %.1f GeV",sqs),nlines,0,nlines,nlines,0,nlines);
243 
244  config_histo2d(htag,"Data Mode","Trigger Line");
245  config_histo2d(htagex,"Data Mode","Trigger Line");
246  config_histo2d(hall);
247 
248 
249  // ********
250  // signal efficiencies
251  // ********
252  TH1F *hsig=new TH1F("hsig"+sqsstr,TString::Format("Signal efficiencies @ %.1f GeV",sqs),nlines,0,nlines);
253  TH1F *hsigi=new TH1F("hsigi"+sqsstr,TString::Format("Signal efficiencies @ %.1f GeV",sqs),nlines,0,nlines);
254  TH1F *hsign=new TH1F("hsign"+sqsstr,TString::Format("Signal efficiencies @ %.1f GeV",sqs),nlines,0,nlines);
255 
256  config_histo1d(hsig, "efficiency [%]", "Data Mode",hsig->GetLineColor(),hsig->GetLineColor(),3003);
257  config_histo1d(hsigi, "efficiency [%]","Data Mode", 2, 2, 3003);
258  config_histo1d(hsign, "efficiency [%]","Data Mode");
259 
260  // ********
261  // background contributions
262  // ********
263  TH1F *hbg=new TH1F("hbg"+sqsstr,TString::Format("Background fractions @ %.1f GeV",sqs),nlines,0,nlines);
264  config_histo1d(hbg, "acc. background [%]", "Trigger Line",1,1,3003,"T");
265  hbg->GetXaxis()->SetBinLabel(nlines,"any");
266 
267  int N = t->GetEntries();
268  int Nbg = 0;
269 
270  c3->cd(); config_pad(gPad); gPad->SetLogz();
271  c1->cd(); config_pad(gPad); gPad->SetLogz();
272 
273  for (int i=0;i<N;++i)
274  {
275  if (fact>1 && i%fact) continue;
276  t->GetEvent(i);
277 
278  int smode = mode;
279  if (smode>1e6) smode/=1000;
280  if (smode>1000) smode%=1000;
281 
282  if (smode==900) Nbg++;
283 
284  hsign->Fill(codeidx[smode]);
285 
286  for (int j=0;j<nlines;++j)
287  hall->Fill(j,nlines-codeidx[smode]-1);
288 
289  if (i%100000==0)
290  {
291  cout <<i<<"/"<<N<<endl;
292 
293  c1->cd(); htag->Draw("colz"); c1->Update();
294  c3->cd(); htagex->Draw("colz"); c3->Update();
295  }
296 
297  if (tag>0)
298  {
299  // simultaneous tag
300  htag->Fill(nlines-1,nlines-codeidx[smode]-1);
301  hsig->Fill(codeidx[smode]);
302  if (smode==900) hbg->Fill(nlines-1);
303 
304  // individual tags
305  for (int j=0;j<nlines-1;++j)
306  {
307  //if (codes[j]==410) continue;
308  if (tagline[j]>0)
309  {
310  htag->Fill(j,nlines-codeidx[smode]-1);
311  if (tagline[j]==tagall) htagex->Fill(j,nlines-codeidx[smode]-1);
312  if (smode==900) hbg->Fill(j);
313  if (j==codeidx[smode]) hsigi->Fill(j);
314  }
315  }
316  }
317  }
318 
319  // ********
320  // matrix of tag line vs data mode
321  // ********
322  c1->cd();
323  htag->Divide(hall); htag->SetMinimum(0.00001); htag->SetMaximum(1.0);
324  htag->Draw("colz");
325 
326  resizePalette(htag);
327  htag->Draw("colz");
328 
329  c1->Update();
330 
331  // ********
332  // matrix of tag line vs data mode, containing the exclusive triggers (= events not tagged by any other line) for certain data modes
333  // ********
334  c3->cd();
335  htagex->Divide(hall); htagex->SetMinimum(0.00001); htagex->SetMaximum(1.0);
336  htagex->Draw("colz");
337 
338  resizePalette(htagex);
339  htagex->Draw("colz");
340 
341  c3->Update();
342 
343  gStyle->SetTitleH(0.06);
344  c2->cd(1);
345  config_pad(gPad,0.14,0.03,0.08,0.08);
346  hsig->Divide(hsign);
347  hsig->Scale(100.);
348  hsig->SetMaximum(100.);
349 
350  hsigi->Divide(hsign);
351  hsigi->Scale(100.);
352  hsigi->SetMaximum(100.);
353 
354  hsig->Draw();
355  hsigi->Draw("same");
356 
357  TLegend *leg=new TLegend(0.08,0.78,0.3,0.92);
358  leg->AddEntry("hsigi"+sqsstr,"dedicated trigger","l");
359  leg->AddEntry("hsig"+sqsstr,"total","l");
360  leg->Draw("same");
361 
362  c2->cd(2);
363  config_pad(gPad,0.14,0.03,0.08,0.08);
364  hbg->Scale(100./Nbg);
365  hbg->Draw();
366 
367  TLine l;
368  l.SetLineColor(1);
369  l.SetLineWidth(2);
370 
371  for (int i=1;i<nlines;++i)
372  {
373  if (drawline.find(codes[i])!=drawline.end() || i==nlines-1)
374  {
375  c1->cd();
376  l.SetLineWidth(2);
377  l.DrawLine(i,0,i,nlines);
378  l.DrawLine(0,nlines-i,nlines,nlines-i);
379  l.SetLineWidth(1);
380  c3->cd();
381  l.SetLineWidth(2);
382  l.DrawLine(i,0,i,nlines);
383  l.DrawLine(0,nlines-i,nlines,nlines-i);
384  l.SetLineWidth(1);
385  c2->cd(1);
386  l.DrawLine(i,0,i,hsig->GetMaximum()*0.75);
387  c2->cd(2);
388  l.DrawLine(i,0,i,hbg->GetMaximum()*0.75);
389  }
390  }
391  c1->Update();
392  c2->Update();
393  c3->Update();
394 
395  if (saveplots>0)
396  {
397  c1->SaveAs(TString::Format("cross_tag_%d.pdf",int(sqs*100)));
398  c1->SaveAs(TString::Format("cross_tag_%d.gif",int(sqs*100)));
399  c2->SaveAs(TString::Format("cross_tag_sigbg_%d.pdf",int(sqs*100)));
400  c2->SaveAs(TString::Format("cross_tag_sigbg_%d.gif",int(sqs*100)));
401 
402  htag->Write();
403  hall->Write();
404  hsig->Write();
405  hsign->Write();
406  hsigi->Write();
407  hbg->Write();
408  ff->Close();
409  }
410  return 0;
411 }
412 
setStyle
void setStyle()
Definition: crosstag.C:39
p
Double_t p
Definition: anasim.C:58
tagm
Int_t tagm
Definition: crosstag.C:23
r
double r
Definition: RiemannTest.C:14
i
Int_t i
Definition: run_full.C:25
m
__m128 m
Definition: P4_F32vec4.h:28
b
TTree * b
Definition: GammaSpectraAnalysis.C:24
tag
Int_t tag
Definition: crosstag.C:23
nlines
Int_t nlines
Definition: crosstag.C:23
crosstag
int crosstag(TString fname, int fact=1, int saveplots=0, double sqs=-1.)
Definition: crosstag.C:203
c2
c2
Definition: plot_dirc.C:39
codes
Int_t codes[60]
Definition: crosstag.C:23
v
__m128 v
Definition: P4_F32vec4.h:4
palette2
void palette2()
Definition: checkphsp2_2.C:56
init
static int init
Definition: ranlxd.cxx:374
mode
Int_t mode
Definition: autocutx.C:47
Double_t
Double_t
Definition: SimCompleteLinkDef.h:6
config_histo
void config_histo(TH1 *h, TString tx, TString ty)
Definition: ana_dsdsj_full.C:19
resizePalette
void resizePalette(TH1 *h)
Definition: crosstag.C:53
f
TFile * f
Definition: bump_analys.C:12
tagall
Int_t tagall
Definition: crosstag.C:23
c1
c1
Definition: plot_dirc.C:35
codeidx
std::map< int, int > codeidx
Definition: crosstag.C:21
c3
c3
Definition: plot_dirc.C:50
config_pad
void config_pad(TVirtualPad *p, double b=0.1, double r=0.12, double t=0.08, double l=0.1)
Definition: checkphsp2_2.C:103
config_histo2d
void config_histo2d(TH1 *h, TString titley="", TString titlex="")
Definition: crosstag.C:159
drawline
std::map< int, int > drawline
Definition: crosstag.C:21
leg
leg
Definition: plot_acceptance_stt.C:24
t
TTree * t
Definition: bump_analys.C:13
config_histo1d
void config_histo1d(TH1 *h, TString titley="", TString titlex="", int lincol=1, int fillcol=0, int fillstyle=0, TString label="M")
Definition: crosstag.C:178
recmode
Int_t recmode
Definition: crosstag.C:23
tagline
Int_t tagline[60]
Definition: crosstag.C:23
Generated on Wed Apr 3 2019 08:43:50 for FairRoot/PandaRoot by   doxygen 1.8.5