don't sort clusters after local reco, do this in AliITSUTrackerGlo
[u/mrichter/AliRoot.git] / ITS / UPGRADE / PrepSummaryKMC.C
CommitLineData
c71b333a 1#if !defined(__CINT__) || defined(__MAKECINT__)
2#include <TROOT.h>
3#include <TFile.h>
4#include <TDirectory.h>
5#include <TF1.h>
6#include <TH1.h>
7#include <TObjArray.h>
8#include <TGraphErrors.h>
9#include "KMCDetector.h"
10#endif
11
12// Process the output file of testDetKMC.C macro to prepare summary
13// for the "summary class" icl (definition of acceptable track).
14// See ProcessSummary comment for details.
15// The names of the graphs are prefixed by `pref`.
16// If outF is provided, the graphs will be stored there
17
18
19TObjArray* PrepSummaryKMC(const char* sumf, int cls, const char* pref=0, const char* outF=0);
20TObjArray* ProcessSummary(TObjArray* sums, int icl, const char* pref=0);
21
22enum {kSigAD,kSigAZ,kSigAPt,kSigD,kSigZ,kSigPt,kEff,kUpd};
23TF1* gs = 0;
24
25TObjArray* PrepSummaryKMC(const char* sumf, int cls, const char* pref, const char* outF)
26{
27 if (!gROOT->GetClass("KMCDetector")) gROOT->LoadMacro("KMCDetector.cxx+");
28 TFile* fl = TFile::Open(sumf);
29 if (!fl) {printf("No file %s\n",sumf); return 0;}
30 TObjArray* arrs = (TObjArray*)gDirectory->Get("trSum");
31 if (!arrs) {printf("No summary in file %s\n",sumf); return 0;}
32 //
33 TObjArray* sums = ProcessSummary(arrs,cls,pref);
34 if (!sums) return 0;
35 //
36 if (outF) {
37 TFile* flOut = TFile::Open(outF,"update");
38 if (!flOut) {printf("Failed to open output file %s\n",outF);}
39 else {
40 flOut->WriteObject(sums,sums->GetName(),"kSingleKey");
41 flOut->Close();
42 delete flOut;
43 printf("Stored array %s in %s\n",sums->GetName(), outF);
44 }
45 }
46 //
47 return sums;
48}
49
50
51TObjArray* ProcessSummary(TObjArray* arrs, int icl, const char* pref)
52{
53 // Process TObjArray (e.g. for set of pt bins) of TObjArray of KMCTrackSummary objects:
54 // pick the KMCTrackSummary for "summary class" icl (definition of acceptable track) and create
55 // graphs vs bin.
56 // These graphs are returned in new TObjArray
57 //
58 TString prefs = pref;
59 if (!gs) gs = new TF1("gs","gaus",-1,1);
60 //
61 int nb = arrs->GetEntriesFast();
62 TObjArray* sums = (TObjArray*) arrs->At(0);
63 int nclass = sums->GetEntriesFast();
64 if (icl>=nclass) {printf("summary set has %d classes only, %d requested\n",nclass,icl);return 0;}
65 //
66 KMCTrackSummary* sm = (KMCTrackSummary*)sums->At(icl);
67 //
68 TH1* h;
69 //
70 h = sm->GetHMCSigDCARPhi(); // MC resolution in transverse DCA
71 TGraphErrors * grSigD = 0;
72 if (h) {
73 grSigD = new TGraphErrors(nb);
74 grSigD->SetName(Form("%s%s",prefs.Data(),h->GetName()));
75 grSigD->SetTitle(Form("%s%s",prefs.Data(),h->GetTitle()));
76 }
77 //
78 TGraphErrors * grSigZ = 0;
79 h = sm->GetHMCSigDCAZ(); // MC resolution in Z DCA
80 if (h) {
81 grSigZ = new TGraphErrors(nb);
82 grSigZ->SetName(Form("%s%s",prefs.Data(),h->GetName()));
83 grSigZ->SetTitle(Form("%s%s",prefs.Data(),h->GetTitle()));
84 }
85 //
86 TGraphErrors * grSigAD = 0; // anaitical estimate for resolution in transverse DCA
87 {
88 grSigAD = new TGraphErrors(nb);
89 grSigAD->SetName(Form("%s%s",prefs.Data(),"sigmaDan"));
90 grSigAD->SetTitle(Form("%s%s",prefs.Data(),"#sigmaD an"));
91 }
92 //
93 TGraphErrors * grSigAZ = 0; // anaitical estimate for resolution in Z DCA
94 {
95 grSigAZ = new TGraphErrors(nb);
96 grSigAZ->SetName(Form("%s%s",prefs.Data(),"sigmaZan"));
97 grSigAZ->SetTitle(Form("%s%s",prefs.Data(),"#sigmaZ an"));
98 }
99 //
100 //
101 TGraphErrors * grSigPt = 0; // MC res. in pt
102 {
103 grSigPt = new TGraphErrors(nb);
104 grSigPt->SetName(Form("%s%s",prefs.Data(),"sigmaPt"));
105 grSigPt->SetTitle(Form("%s%s",prefs.Data(),"#sigmaPt"));
106 }
107 //
108 TGraphErrors * grSigAPt = 0; // analitycal res. in pt
109 {
110 grSigAPt = new TGraphErrors(nb);
111 grSigAPt->SetName(Form("%s%s",prefs.Data(),"sigmaPtan"));
112 grSigAPt->SetTitle(Form("%s%s",prefs.Data(),"#sigmaPt an"));
113 }
114
115 //
116 TGraphErrors * grEff = 0; // MC efficiency
117 {
118 grEff = new TGraphErrors(nb);
119 grEff->SetName(Form("%s_rate",prefs.Data()));
120 grEff->SetTitle(Form("%s Rate",prefs.Data()));
121 }
122 //
123 TGraphErrors * grUpd = 0; // number of Kalman track updates
124 {
125 grUpd = new TGraphErrors(nb);
126 grUpd->SetName(Form("%s_updCalls",prefs.Data()));
127 grUpd->SetTitle(Form("%s Updates",prefs.Data()));
128 }
129 //
130 for (int ib=0;ib<nb;ib++) {
131 sums = (TObjArray*) arrs->At(ib);
132 sm = (KMCTrackSummary*)sums->At(icl);
133 KMCProbe& prbRef = sm->GetRefProbe();
134 KMCProbe& prbAn = sm->GetAnProbe();
135
136 double pt = prbRef.Pt();
137 //
138 if (grSigAD) {
139 grSigAD->SetPoint(ib, pt,prbAn.GetSigmaY2()>0 ? TMath::Sqrt(prbAn.GetSigmaY2()) : 0.);
140 }
141 //
142 if (grSigAZ) {
143 grSigAZ->SetPoint(ib, pt,prbAn.GetSigmaZ2()>0 ? TMath::Sqrt(prbAn.GetSigmaZ2()) : 0.);
144 }
145 //
146 if (grSigAPt) {
147 double pts = TMath::Sqrt(prbAn.GetSigma1Pt2());
148 grSigAPt->SetPoint(ib, pt,pts>0 ? pts*pt : 0.);
149 }
150 //
151 if (grSigPt) {
152 h = sm->GetHMCSigPt();
153 h->Fit(gs,"0q");
154 grSigPt->SetPoint(ib, pt, gs->GetParameter(2));
155 grSigPt->SetPointError(ib, 0, gs->GetParError(2));
156 }
157 //
158 if (grSigD) {
159 h = sm->GetHMCSigDCARPhi();
160 h->Fit(gs,"0q");
161 grSigD->SetPoint(ib, pt,gs->GetParameter(2));
162 grSigD->SetPointError(ib, 0,gs->GetParError(2));
163 }
164 //
165 if (grSigZ) {
166 h = sm->GetHMCSigDCAZ();
167 h->Fit(gs,"0q");
168 grSigZ->SetPoint(ib, pt,gs->GetParameter(2));
169 grSigZ->SetPointError(ib, 0,gs->GetParError(2));
170 }
171 //
172 if (grEff) {
173 grEff->SetPoint(ib, pt,sm->GetEff());
174 grEff->SetPointError(ib, 0,sm->GetEffErr());
175 }
176 //
177 if (grUpd) {
178 grUpd->SetPoint(ib, pt,sm->GetUpdCalls());
179 grUpd->SetPointError(ib, 0, 0);
180 }
181 }
182 //
183 TObjArray* dest = new TObjArray();
184 dest->AddAtAndExpand(grSigAD,kSigAD);
185 dest->AddAtAndExpand(grSigAZ,kSigAZ);
186 dest->AddAtAndExpand(grSigAPt,kSigAPt);
187 dest->AddAtAndExpand(grSigD,kSigD);
188 dest->AddAtAndExpand(grSigZ,kSigZ);
189 dest->AddAtAndExpand(grSigPt,kSigPt);
190 dest->AddAtAndExpand(grEff,kEff);
191 dest->AddAtAndExpand(grUpd,kUpd);
192 //
193 if (!prefs.IsNull()) dest->SetName(pref);
194 return dest;
195}