]>
Commit | Line | Data |
---|---|---|
9dc17dc7 | 1 | /* |
2 | ||
3 | //load mag field | |
4 | AliRunLoader *loader = AliRunLoader::Open("galice.root"); | |
5 | loader->LoadgAlice(); | |
6 | gAlice = loader->GetAliRun(); | |
7 | AliMagF * field = gAlice->Field(); | |
e2cd3644 | 8 | AliTracker::SetFieldMap(gAlice->Field(),1); |
9dc17dc7 | 9 | TFile fg("geometry.root"); |
10 | gGeoManager = (TGeoManager*)fg.Get("Geo"); | |
7dd16671 | 11 | if (!gGeoManager) TGeoManager::Import("geometry.root"); |
9dc17dc7 | 12 | .L AliGenInfo.C+ |
13 | .L AliESDComparisonMI.C+ | |
7dd16671 | 14 | .L AliTRDComparison.C+g |
9dc17dc7 | 15 | //.L AliESDComparisonPIC.C+ |
16 | .x TDR_style.C | |
17 | ||
18 | MakeCompTr(); | |
19 | MakeTree(); | |
20 | ||
e2cd3644 | 21 | |
7dd16671 | 22 | |
9dc17dc7 | 23 | MakeComp(); |
24 | MakeAlias(); | |
25 | TFile fff("TRDdebug.root"); | |
26 | TTree * tree =(TTree*)fff.Get("tracklet"); | |
27 | TTree * tree2 =(TTree*)fff.Get("Tracks"); | |
7dd16671 | 28 | TTree * treefind =(TTree*)fff.Get("Find"); |
29 | TTree * treebug =(TTree*)fff.Get("Bug"); | |
30 | TTree * treetofs =(TTree*)fff.Get("tofseed"); | |
31 | AliComparisonDraw compfind; | |
32 | compfind.fTree = treefind; | |
33 | AliComparisonDraw comptr; | |
34 | comptr.fTree = tree; | |
9dc17dc7 | 35 | |
7dd16671 | 36 | TFile fff2("TOFdebug.root"); |
37 | TTree * treetof =(TTree*)fff2.Get("Info"); | |
38 | TTree * treetoftrd =(TTree*)fff2.Get("TRD"); | |
9dc17dc7 | 39 | |
7dd16671 | 40 | AliComparisonDraw comp2; |
41 | comp2->fTree =tree; | |
42 | AliComparisonDraw comptof; | |
43 | comptof->fTree =treetof; | |
9dc17dc7 | 44 | |
45 | ||
46 | */ | |
47 | ||
48 | //ROOT includes | |
49 | #include "Rtypes.h" | |
50 | #include "TFile.h" | |
51 | #include "TTree.h" | |
52 | #include "TChain.h" | |
53 | #include "TCut.h" | |
54 | #include "TString.h" | |
55 | #include "TBenchmark.h" | |
56 | #include "TStopwatch.h" | |
57 | #include "TParticle.h" | |
58 | #include "TSystem.h" | |
59 | #include "TTimer.h" | |
60 | #include "TVector3.h" | |
61 | #include "TH1F.h" | |
62 | #include "TH2F.h" | |
63 | #include "TCanvas.h" | |
64 | #include "TPad.h" | |
65 | #include "TF1.h" | |
66 | #include "TView.h" | |
67 | #include "TGeometry.h" | |
68 | #include "TPolyLine3D.h" | |
69 | #include "TLegend.h" | |
7dd16671 | 70 | #include "TGraph.h" |
9dc17dc7 | 71 | //ALIROOT includes |
72 | #include "AliRun.h" | |
73 | #include "AliStack.h" | |
74 | #include "AliSimDigits.h" | |
75 | #include "AliTPCParam.h" | |
76 | #include "AliTPC.h" | |
77 | #include "AliTPCLoader.h" | |
78 | #include "AliDetector.h" | |
79 | #include "AliTrackReference.h" | |
80 | #include "AliTPCParamSR.h" | |
81 | #include "AliTracker.h" | |
82 | #include "AliMagF.h" | |
83 | #include "AliHelix.h" | |
84 | #include "AliPoints.h" | |
85 | #include "AliESDtrack.h" | |
86 | #include "AliTRDtrack.h" | |
87 | #include "AliITStrackMI.h" | |
88 | #include "AliESDV0MI.h" | |
89 | #include "AliESDkink.h" | |
7014e413 | 90 | //#include "AliTRDparameter.h" |
9dc17dc7 | 91 | #include "AliTRDgeometryFull.h" |
92 | #include "AliTRDcluster.h" | |
7dd16671 | 93 | #include "AliTRDtracker.h" |
94 | #include "AliTOFtrack.h" | |
95 | #include "../TOF/AliTOFtrackerMI.h" | |
9dc17dc7 | 96 | // COMPARISON includes |
97 | #include "../STEER/AliGenInfo.h" | |
98 | #include "../STEER/AliESDComparisonMI.h" | |
99 | ||
100 | ||
101 | ||
102 | AliComparisonDraw comp; | |
103 | AliComparisonDraw comptr; | |
7014e413 | 104 | //AliTRDparameter par("pica","vyjebana"); |
9dc17dc7 | 105 | AliTRDgeometryFull geom; |
106 | AliRunLoader *fLoader=0; //= AliRunLoader::Open(fnGalice); | |
107 | ||
7dd16671 | 108 | TCut cnesd("cnesd","abs(RC.fLabel)!=MC.fLabel"); |
9dc17dc7 | 109 | TCut cprim("cprim","TMath::Sqrt(MC.fVDist[0]**2+MC.fVDist[1]**2)<0.01&&abs(MC.fVDist[2])<0.01"); |
110 | TCut citsin("cprim","TMath::Sqrt(MC.fVDist[0]**2+MC.fVDist[1]**2)<0.01&&abs(MC.fVDist[2])<5"); | |
111 | TCut clength("clength","abs(TRD0.fIntegratedLength-MC.fTOFReferences[0].fLength)<30"); | |
112 | TCut c1("c1","TRD1.fStopped==0"); // track propagated to refernce | |
113 | TCut c0("c0","TRD0.fStopped==0"); // track propagated to reference | |
114 | TCut c2("c2","MC.fNTOFRef>0"); //track in tof | |
115 | TCut c3("c3","MC.fNTOFRef>0&&MC.fNTRDRef>5"); //track in tof and TRD | |
116 | TH1F hpull("hpull","hpull",50,-10,10); | |
117 | TH1F hpull6("hpull6","hpull6",50,-6,6); | |
118 | ||
119 | ||
120 | class AliTRDInfoMC: public TObject{ | |
121 | public: | |
122 | void Update(AliMCInfo *info); // update TRD info | |
123 | void Reset(); // clear TRD info | |
124 | Short_t fLayer[6]; //indicates layer is on | |
125 | Short_t fLayerA[6]; // indicates layer after is ON | |
126 | Short_t fNLayers; // number of crossed layers | |
127 | Short_t fLast; // last crossed layer | |
128 | AliTrackReference fRef[6]; // trd reference | |
7dd16671 | 129 | AliTrackReference fRefLocal[6]; // trd local reference |
9dc17dc7 | 130 | ClassDef(AliTRDInfoMC,1) // TRD MC information |
131 | }; | |
132 | ||
9dc17dc7 | 133 | |
134 | class AliTRDInfoCl: public TObject{ | |
135 | public: | |
7dd16671 | 136 | AliTRDInfoCl(); |
137 | void Update(); | |
138 | AliTRDcluster fCl[6][25]; // clusters | |
139 | Float_t fCount[6][25]; // number of clusters of given label | |
140 | Int_t fNClusters; // number of clusters | |
141 | Float_t fMeanQ[6]; // mean Q | |
142 | Int_t fMaxBin[6]; // max time bin | |
143 | Int_t fNcl[6]; // number of clusters | |
144 | Float_t fMeanQS[6]; // mean Q | |
145 | Int_t fMaxBinS[6]; // max time bin | |
146 | ClassDef(AliTRDInfoCl,1) // TRD MC information | |
9dc17dc7 | 147 | }; |
148 | ||
149 | ClassImp(AliTRDInfoCl) | |
150 | ||
151 | ||
152 | ||
7dd16671 | 153 | class AliTRDInfoSeed: public TObject{ |
154 | public: | |
155 | AliTRDInfoSeed(){fLabel=-1; fLabel1=-1;fLabel2=-1;} | |
156 | void Update(AliTRDseed *seeds); | |
157 | AliTRDseed fSeeds[6]; // seeds from seeding | |
158 | Float_t fFakeRatio; // | |
159 | Int_t fLabel; // seed label | |
160 | Int_t fLabel1; // full track label | |
161 | Int_t fLabel2; // second full track label | |
162 | Int_t fNLayers; // | |
163 | Int_t fNUsed; // | |
164 | ClassDef(AliTRDInfoSeed,1) // TRD MC information | |
165 | }; | |
166 | ||
167 | ||
168 | ClassImp(AliTRDInfoMC) | |
169 | ClassImp(AliTRDInfoCl) | |
170 | ClassImp(AliTRDInfoSeed) | |
171 | ||
172 | void AliTRDInfoSeed::Update(AliTRDseed *seeds){ | |
173 | // | |
174 | // | |
175 | // | |
176 | fNLayers=0; | |
177 | for (Int_t i=0;i<6;i++){ | |
178 | fSeeds[i] = seeds[i]; | |
179 | if (seeds[i].fN2>10) fNLayers++; | |
180 | } | |
181 | } | |
182 | ||
9dc17dc7 | 183 | void AliTRDInfoMC::Reset(){ |
184 | // | |
185 | // | |
186 | fNLayers=0; | |
187 | fLast =0; | |
188 | for (Int_t ilayer = 0; ilayer<6; ilayer++){ | |
189 | fLayer[ilayer]=0; | |
190 | fLayerA[ilayer]=0; | |
191 | } | |
192 | } | |
193 | void AliTRDInfoMC::Update(AliMCInfo *info) | |
194 | { | |
195 | // | |
196 | // | |
197 | Reset(); | |
9dc17dc7 | 198 | for (Int_t i=0;i<info->fNTRDRef;i++){ |
199 | AliTrackReference * ref = (AliTrackReference*)info->fTRDReferences->At(i); | |
7dd16671 | 200 | Float_t localX = ref->LocalX(); |
9dc17dc7 | 201 | for (Int_t jplane = 0; jplane<6; jplane++){ |
7dd16671 | 202 | if (TMath::Abs(localX - (300+jplane*14.))<6){ |
9dc17dc7 | 203 | fLayer[jplane]++; |
204 | for (Int_t kplane=jplane;kplane>=0;kplane--){ | |
205 | fLayerA[kplane]++; | |
206 | fRef[jplane]=*ref; | |
7dd16671 | 207 | // |
208 | fRefLocal[jplane].SetTrack(ref->GetTrack()); | |
209 | fRefLocal[jplane].SetPosition(ref->LocalX(), ref->LocalY(),ref->Z()); | |
210 | Float_t alpha = ref->Alpha(); | |
211 | Float_t dir[3] = {ref->Px(), ref->Py(), ref->Pz()}; | |
212 | Float_t dirR[3]; | |
213 | dirR[0] = dir[0]*TMath::Cos(-alpha) - dir[1]*TMath::Sin(-alpha); | |
214 | dirR[1] = dir[0]*TMath::Sin(-alpha) + dir[1]*TMath::Cos(-alpha); | |
215 | fRefLocal[jplane].SetMomentum(1, dirR[1]/dirR[0],dir[2]/dir[0]); | |
9dc17dc7 | 216 | } |
217 | } | |
218 | } | |
219 | } | |
220 | for (Int_t jplane = 0; jplane<6; jplane++){ | |
221 | if (fLayer[jplane]>0) fNLayers++; | |
222 | if (fLayer[jplane]>0) fLast = jplane; | |
223 | } | |
224 | } | |
225 | ||
7dd16671 | 226 | AliTRDInfoCl::AliTRDInfoCl(){ |
227 | // | |
228 | // | |
229 | fNClusters =0; | |
230 | for (Int_t ilayer=0;ilayer<6;ilayer++){ | |
231 | fMeanQ[ilayer]=0; | |
232 | fMaxBin[ilayer]=0; | |
233 | fNcl[ilayer]=0; | |
234 | fMeanQS[ilayer]=0; | |
235 | fMaxBinS[ilayer]=0; | |
236 | for (Int_t itime=0;itime<25;itime++){ | |
237 | fCount[ilayer][itime]=0; | |
238 | } | |
239 | } | |
240 | } | |
241 | ||
242 | ||
243 | void AliTRDInfoCl::Update(){ | |
244 | fNClusters = 0; | |
245 | for (Int_t iplane=0;iplane<6;iplane++){ | |
246 | fMeanQ[iplane]=-1; | |
247 | fNcl[iplane] = 0; | |
248 | fMaxBin[iplane]=-1; | |
249 | // | |
250 | fMeanQS[iplane]=-1; | |
251 | fMaxBinS[iplane]=-1; | |
252 | ||
253 | Float_t maxQ =-10; | |
254 | Float_t maxQS =-10; | |
255 | Int_t sums =0; | |
256 | for (Int_t itime=0;itime<25;itime++){ | |
257 | if (fCl[iplane][itime].GetQ()>2){ | |
258 | fNcl[iplane]++; | |
259 | fMeanQ[iplane]+=fCl[iplane][itime].GetQ(); | |
260 | if (fCl[iplane][itime].GetQ()>=maxQ){ | |
261 | fMaxBin[iplane]= itime; | |
262 | maxQ = fCl[iplane][itime].GetQ(); | |
263 | } | |
264 | } | |
265 | if (fCl[iplane][itime].GetSumS()>2){ | |
266 | sums++; | |
267 | fMeanQS[iplane]+=fCl[iplane][itime].GetSumS(); | |
268 | if (fCl[iplane][itime].GetSumS()>=maxQS){ | |
269 | fMaxBinS[iplane]= itime; | |
270 | maxQS = fCl[iplane][itime].GetSumS(); | |
271 | } | |
272 | } | |
273 | } | |
274 | if (fNcl[iplane]>0) fMeanQ[iplane]/=Float_t(fNcl[iplane]); | |
275 | if (sums>0) fMeanQS[iplane]/=Float_t(sums); | |
276 | fNClusters+=fNcl[iplane]; | |
277 | } | |
278 | } | |
279 | ||
9dc17dc7 | 280 | |
281 | ||
282 | ||
283 | void MaxBinNorm(const char * expr = "TRD0.fTimBinPlane"){ | |
284 | // | |
285 | // max time bin probability distribution | |
286 | // | |
287 | char var[100]; | |
288 | sprintf(var,"%s>>hispi",expr); | |
289 | TH1F * hispi= new TH1F("hispi","hispi",20,1,22); | |
290 | comp.fTree->Draw(var,"abs(MC.fPdg)==211&&TRD0.fN>90&&MC.fParticle.P()>3"); | |
291 | hispi->Scale(1./(hispi->GetEntries())); | |
292 | hispi->Draw(); | |
293 | sprintf(var,"%s>>hisel",expr); | |
294 | TH1F *hisel = new TH1F("hisel","hisel",20,1,22); | |
295 | comp.fTree->Draw(var,"abs(MC.fPdg)==11&&TRD0.fN>90&&MC.fParticle.P()>3"); | |
296 | hisel->Scale(1./(hisel->GetEntries())); | |
297 | hispi->SetMarkerStyle(25); hispi->SetLineStyle(2); | |
298 | hisel->SetMarkerStyle(26); hisel->SetLineStyle(1); | |
299 | hisel->SetXTitle("Maximal Time Bin"); | |
300 | hisel->SetYTitle("Probability"); | |
301 | hisel->Draw(); | |
302 | hispi->Draw("same"); | |
303 | TLegend *legend = new TLegend(0.55,0.12,0.85,0.35, "Maximal Time Bin position"); | |
304 | legend->SetBorderSize(1); | |
305 | legend->AddEntry(hispi ,"Pions"); | |
306 | legend->AddEntry(hisel ,"Electrons"); | |
307 | legend->Draw(); | |
308 | ||
309 | } | |
310 | ||
311 | ||
312 | void dedxNorm(const char * expr = "TRD0.fdEdxPlane/20"){ | |
313 | // | |
314 | // max time bin probability distribution | |
315 | // | |
316 | char var[100]; | |
317 | sprintf(var,"%s>>hispi",expr); | |
318 | TH1F * hispi= new TH1F("hispi","hispi",100,1,120); | |
319 | comp.fTree->Draw(var,"abs(MC.fPdg)==211&&TRD0.fN>90&&MC.fParticle.P()>3"); | |
320 | hispi->Scale(1./(hispi->GetEntries())); | |
321 | hispi->Draw(); | |
322 | sprintf(var,"%s>>hisel",expr); | |
323 | TH1F *hisel = new TH1F("hisel","hisel",100,1,120); | |
324 | comp.fTree->Draw(var,"abs(MC.fPdg)==11&&TRD0.fN>90&&MC.fParticle.P()>3"); | |
325 | hisel->Scale(1./(hisel->GetEntries())); | |
326 | hispi->SetMarkerStyle(25); hispi->SetLineStyle(2); | |
327 | hisel->SetMarkerStyle(26); hisel->SetLineStyle(1); | |
328 | hispi->SetXTitle("Normalized dEdx"); | |
329 | hispi->SetYTitle("Probability"); | |
330 | hispi->Draw(); | |
331 | hisel->Draw("same"); | |
332 | TLegend *legend = new TLegend(0.55,0.12,0.85,0.35, "Normalized dEdx"); | |
333 | legend->SetBorderSize(1); | |
334 | legend->AddEntry(hispi ,"Pions"); | |
335 | legend->AddEntry(hisel ,"Electrons"); | |
336 | legend->Draw(); | |
337 | } | |
338 | ||
339 | ||
340 | void EfficiencyPt() | |
341 | { | |
342 | // | |
343 | // | |
344 | // | |
345 | comp.Eff("TR.Pt()","InfoMC.fNLayers>3&&abs(MC.fPdg)!=11"+c2+cprim,"TRD0.fN/InfoMC.fNLayers>10&&TRD0.fLab>0",10,0.3,2); | |
346 | TH1F * his = (TH1F*)comp.fRes->Clone(); | |
347 | comp.Eff("TR.Pt()","InfoMC.fNLayers>3&&abs(MC.fPdg)!=11"+c2+cprim,"TRD0.fN/InfoMC.fNLayers>10&&TRD0.fLab<0",10,0.3,2); | |
348 | TH1F * hisf = (TH1F*)comp.fRes->Clone(); | |
349 | his->SetMarkerStyle(25); his->SetLineStyle(2); | |
350 | hisf->SetMarkerStyle(26); hisf->SetLineStyle(1); | |
351 | his->SetXTitle("P_{t}(GeV/c)"); | |
352 | his->SetYTitle("Tracking efficiency (%)"); | |
353 | his->Draw(); | |
354 | hisf->Draw("same"); | |
355 | TLegend *legend = new TLegend(0.55,0.16,0.85,0.4, "TRD efficciency"); | |
356 | legend->SetBorderSize(1); | |
357 | legend->AddEntry(his ,"Efficiency"); | |
358 | legend->AddEntry(hisf ,"Fake tracks"); | |
359 | legend->Draw(); | |
360 | ||
361 | } | |
362 | ||
363 | void EfficiencyDip() | |
364 | { | |
365 | // | |
366 | // | |
367 | // | |
368 | comp.Eff("atan(TR.fPz/TR.Pt())","InfoMC.fNLayers>3&&abs(MC.fPdg)!=11&&TR.Pt()>0.4"+c2+cprim,"TRD0.fN/InfoMC.fNLayers>10&&TRD0.fLab>0",10,-0.9,0.9); | |
369 | TH1F * his = (TH1F*)comp.fRes->Clone(); | |
370 | comp.Eff("atan(TR.fPz/TR.Pt())","InfoMC.fNLayers>3&&abs(MC.fPdg)!=11&&TR.Pt()>0.4"+c2+cprim,"TRD0.fN/InfoMC.fNLayers>10&&TRD0.fLab<0",10,-0.9,0.9); | |
371 | TH1F * hisf = (TH1F*)comp.fRes->Clone(); | |
372 | his->SetMarkerStyle(25); his->SetLineStyle(2); | |
373 | hisf->SetMarkerStyle(26); hisf->SetLineStyle(1); | |
374 | his->SetXTitle("#lambda(rad)"); | |
375 | his->SetYTitle("Tracking efficiency (%)"); | |
376 | his->Draw(); | |
377 | hisf->Draw("same"); | |
378 | TLegend *legend = new TLegend(0.55,0.16,0.85,0.4, "TRD efficciency"); | |
379 | legend->SetBorderSize(1); | |
380 | legend->AddEntry(his ,"Efficiency"); | |
381 | legend->AddEntry(hisf ,"Fake tracks"); | |
382 | legend->Draw(); | |
383 | } | |
384 | ||
385 | void ClusterRatio() | |
386 | { | |
387 | TH1F *hisr = new TH1F("hisr","Percentage of found clusters",80,0,120); | |
388 | comp.fTree->Draw("100*TRD0.fN/(20*InfoMC.fNLayers)>>hisr","InfoMC.fNLayers>3&&TRD0.fLab>0&&abs(MC.fPdg)!=11&&TR.Pt()>0.4"+cprim+c2); | |
389 | hisr->SetXTitle("(%)"); | |
390 | hisr->Draw(); | |
391 | ||
392 | } | |
393 | ||
7dd16671 | 394 | void DEestimate(){ |
395 | TH1F * hderel = new TH1F("hderel","hderel",100,-150,150); | |
396 | comp.fTree->Draw("100*derel>>hderel",""+c0+c2+cprim); | |
397 | hderel->Fit("gaus"); | |
398 | hderel->SetXTitle("Energy loss correction resolution (%)"); | |
399 | hderel->Draw(); | |
400 | } | |
401 | ||
402 | void DEestimate2(){ | |
403 | // | |
404 | TH1F * hderel0 = new TH1F("hderel0","hderel0",100,-0,30); | |
405 | TH1F * hderel1 = new TH1F("hderel1","hderel1",100,-0,30); | |
406 | TH1F * hderela = new TH1F("hderela","hderela",100,-0,30); | |
407 | comp.fTree->Draw("100*(TPCE-TOFE)/TPCE>>hderela","Pt<1"+c0+c2+cprim); | |
408 | comp.fTree->Draw("100*(TPCE-TOFE)/TPCE>>hderel0","Pt<1&&TRD0.fBudget[0]<10"+c0+c2+cprim); | |
409 | comp.fTree->Draw("100*(TPCE-TOFE)/TPCE>>hderel1","Pt<1&&TRD0.fBudget[0]>10"+c0+c2+cprim); | |
410 | hderela->SetXTitle("(E_{TPC}-E_{TOF})/E_{TPC}(%)"); | |
411 | hderela->SetLineStyle(1); | |
412 | hderel0->SetLineStyle(2); | |
413 | hderel1->SetLineStyle(3); | |
414 | hderela->Draw(""); | |
415 | hderel0->Draw("same"); | |
416 | hderel1->Draw("same"); | |
417 | TLegend *legend = new TLegend(0.55,0.12,0.85,0.35, "Relative Energy Loss Between TPC and TOF"); | |
418 | legend->SetBorderSize(1); | |
419 | legend->AddEntry(hderela ,"All tracks"); | |
420 | legend->AddEntry(hderela ,"Non crossing of boundary"); | |
421 | legend->AddEntry(hderela ,"Crossing boundary"); | |
422 | legend->Draw(); | |
423 | } | |
424 | ||
425 | ||
426 | void EffBudget(){ | |
427 | // | |
428 | // | |
429 | // | |
430 | TF1 f1("f1","100*exp(-[0]*x)",0,100); | |
431 | comp.Eff("TRD.fBudget[0]","abs(ThetaM)<1&&abs(MC.fPdg)==211"+cprim,c2,10,1,100); | |
432 | TH1F * heffb = (TH1F*)comp.fRes->Clone(); | |
433 | heffb->SetXTitle("Material budget (g/cm^{2})"); | |
434 | heffb->SetYTitle("Efficiency to reach TOF (%)"); | |
435 | heffb->Fit(&f1); | |
436 | heffb->Draw(); | |
437 | ||
438 | ||
439 | } | |
440 | ||
441 | void TOFbackround(){ | |
442 | TH1F * his = new TH1F("his","his",100,-100,10000); | |
443 | TCut cut = "InfoTOF.fNCluster>1&&tofm0"; | |
444 | comp.fTree->Draw("InfoTOF.fTime1-33*sqrt(InfoTOF.fCl1.fZ**2+InfoTOF.fCl1.fR**2)>>his","InfoTOF.fTime1<100000"+cut); | |
445 | his->SetXTitle("TOF Time (ns)"); | |
446 | } | |
447 | ||
448 | void TOFsignal(){ | |
449 | TCut cut = "InfoTOF.fNCluster>0&&abs(MC.fPdg)==2212&&P<2"; | |
450 | // | |
451 | TH1F * his0 = new TH1F("his0","his0",100,-1000,1000); | |
452 | comp.fTree->Draw("InfoTOF.fTime0-InfoTOF.fTimes0[4]>>his0","tofm0"+cut); | |
453 | TH1F * hisf = new TH1F("hisf","hisf",100,-1000,1000); | |
454 | comp.fTree->Draw("InfoTOF.fTime0-InfoTOF.fTimes0[4]>>hisf","!(tofm0||tofm0s)"+cut); | |
455 | his0->Draw(); | |
456 | hisf->Draw("same"); | |
457 | his0->SetXTitle("TOF Time (ns)"); | |
458 | } | |
459 | ||
460 | ||
461 | ||
462 | ||
463 | void EffPIDK(Float_t xmin=0.4, Float_t xmax=2){ | |
464 | // | |
465 | // | |
466 | //Float_t xmin=0.4; Float_t xmax=2; | |
467 | // TCut cut= "1"; | |
468 | // // TCut cut = c2; | |
469 | // // TCut cut = "tofm" | |
470 | // AliTOFtrackInfo::SetPIDMatchCuts(); | |
471 | // comp.Eff("P","TRD0.fIntegratedLength>360&&abs(MC.fPdg)==321"+cut,"InfoTOF.GetCombPID(-1,0.1,0)==3",6,xmin,xmax); | |
472 | // TH1F *his0 = (TH1F*)comp.fRes->Clone(); | |
473 | // comp.Eff("P","TRD0.fIntegratedLength>360&&abs(MC.fPdg)==321"+cut,"InfoTOF.GetCombPID(-1,0.1,1)==3",6,xmin,xmax); | |
474 | // TH1F *his1 = (TH1F*)comp.fRes->Clone(); | |
475 | // //AliTOFtrackInfo::SetPIDMatchCuts(-1,-1,-1,100); // disable pid match cuts | |
476 | // //comp.Eff("P","TRD0.fIntegratedLength>360&&abs(MC.fPdg)==321"+cut,"InfoTOF.GetCombPID(-1,0.0,0)==3",6,xmin,xmax); | |
477 | // //TH1F *hisno = (TH1F*)comp.fRes->Clone(); | |
478 | // //AliTOFtrackInfo::SetPIDMatchCuts(); //enable again | |
479 | // // | |
480 | // comp.Eff("P","TRD0.fIntegratedLength>360&&abs(MC.fPdg)==321"+cut,"InfoTOF.GetCombPID(-1,0.1,2)==3",6,xmin,xmax); | |
481 | // TH1F *hisold = (TH1F*)comp.fRes->Clone(); | |
482 | // comp.Eff("P","TRD0.fIntegratedLength>360&&abs(MC.fPdg)==321"+cut,"InfoTOF.GetCombPID(-1,0.0,-1)==3",6,xmin,xmax); | |
483 | // TH1F *histpc = (TH1F*)comp.fRes->Clone(); | |
484 | // his0->SetMarkerStyle(21); | |
485 | // his1->SetMarkerStyle(24); | |
486 | // // hisno->SetMarkerStyle(26); | |
487 | // hisold->SetMarkerStyle(25); | |
488 | // histpc->SetMarkerStyle(23); | |
489 | // // | |
490 | // his0->Draw(); | |
491 | // his1->Draw("same"); | |
492 | // //hisno->Draw("same"); | |
493 | // hisold->Draw("same"); | |
494 | // histpc->Draw("same"); | |
495 | // TLegend *legend = new TLegend(0.55,0.12,0.85,0.35, "PID efficiency "); | |
496 | // legend->SetBorderSize(1); | |
497 | // legend->AddEntry(his0 ,"TPC+New TOF v0"); | |
498 | // legend->AddEntry(his1 ,"TPC+New TOF v1"); | |
499 | // //legend->AddEntry(hisno ,"TPC+New TOF v0 - no PID match cuts"); | |
500 | // legend->AddEntry(hisold ,"TPC+Old TOF"); | |
501 | // legend->AddEntry(histpc ,"TPC only"); | |
502 | // legend->Draw(); | |
503 | // | |
504 | } | |
505 | ||
506 | void EffPIDP(Float_t xmin=0.4, Float_t xmax=2){ | |
507 | // | |
508 | // | |
509 | //Float_t xmin=0.4; Float_t xmax=4; | |
510 | TCut cut= "1"; | |
511 | // TCut cut= c2; | |
512 | //TCut cut = "tofm0" | |
513 | //AliTOFtrackInfo::SetPIDMatchCuts(-1,-1,-1,1000); // default cuts | |
514 | // AliTOFtrackInfo::SetPIDMatchCuts(); // default cuts | |
515 | // comp.Eff("P","TRD0.fIntegratedLength>360&&abs(MC.fPdg)==2212"+cut,"InfoTOF.GetCombPID(-1,0.1,0)==4",6,xmin,xmax); | |
516 | // TH1F *his0 = (TH1F*)comp.fRes->Clone(); | |
517 | // comp.Eff("P","TRD0.fIntegratedLength>360&&abs(MC.fPdg)==2212"+cut,"InfoTOF.GetCombPID(-1,0.1,1)==4",6,xmin,xmax); | |
518 | // TH1F *his1 = (TH1F*)comp.fRes->Clone(); | |
519 | // comp.Eff("P","TRD0.fIntegratedLength>360&&abs(MC.fPdg)==2212"+cut,"InfoTOF.GetCombPID(-1,0.1,2)==4",6,xmin,xmax); | |
520 | // TH1F *hisold = (TH1F*)comp.fRes->Clone(); | |
521 | // comp.Eff("P","TRD0.fIntegratedLength>360&&abs(MC.fPdg)==2212"+cut,"InfoTOF.GetCombPID(-1,0.1,-1)==4",6,xmin,xmax); | |
522 | // TH1F *histpc = (TH1F*)comp.fRes->Clone(); | |
523 | // his0->SetMarkerStyle(21); | |
524 | // his1->SetMarkerStyle(24); | |
525 | // hisold->SetMarkerStyle(25); | |
526 | // histpc->SetMarkerStyle(23); | |
527 | // his0->Draw(); | |
528 | // his1->Draw("same"); | |
529 | // hisold->Draw("same"); | |
530 | // histpc->Draw("same"); | |
531 | // TLegend *legend = new TLegend(0.55,0.12,0.85,0.35, "PID efficiency "); | |
532 | // legend->SetBorderSize(1); | |
533 | // legend->AddEntry(his0 ,"TPC+New TOF v0"); | |
534 | // legend->AddEntry(his1 ,"TPC+New TOF v1"); | |
535 | // legend->AddEntry(hisold ,"TPC+Old TOF"); | |
536 | // legend->AddEntry(histpc ,"TPC only"); | |
537 | // legend->Draw(); | |
538 | ||
539 | // | |
540 | } | |
541 | ||
542 | ||
543 | ||
9dc17dc7 | 544 | |
545 | ||
546 | void PtRes(){ | |
547 | ||
548 | ||
549 | } | |
550 | ||
551 | ||
7dd16671 | 552 | TH1F * MakeCumul(TH1F *his,Bool_t norm = kTRUE) |
9dc17dc7 | 553 | { |
554 | TH1F *hcumul = (TH1F*)his->Clone(); | |
7dd16671 | 555 | char name[1000]; |
556 | sprintf(name,"N%f",gRandom->Rndm()); | |
557 | hcumul->SetName(name); | |
558 | hcumul->SetTitle(name); | |
559 | // | |
9dc17dc7 | 560 | Float_t all = (Float_t)his->GetSum(); |
7dd16671 | 561 | if (!norm) all =1; |
9dc17dc7 | 562 | for (Int_t i=0;i<=his->GetNbinsX();i++){ |
563 | hcumul->SetBinContent(i,his->Integral(0,i)/all); | |
564 | } | |
565 | hcumul->SetFillColor(0); | |
566 | return hcumul; | |
567 | } | |
568 | ||
7dd16671 | 569 | TGraph* His01(TH1F *his0,TH1F *his1){ |
570 | Int_t nbins = his0->GetNbinsX(); | |
571 | Double_t *x = new Double_t[nbins]; | |
572 | Double_t *y = new Double_t[nbins]; | |
573 | for (Int_t i=0;i<nbins;i++) { | |
574 | x[i] = his0->GetBinContent(i+1); | |
575 | y[i] = his1->GetBinContent(i+1); | |
576 | } | |
577 | return new TGraph(nbins,x,y); | |
578 | } | |
579 | ||
580 | ||
9dc17dc7 | 581 | TH1F * PullCumul(const char *var, TCut cut, Int_t div=200, Float_t max=20){ |
582 | TH1F * hpullb = new TH1F("hpullb","hpullb",div,0,max); | |
583 | char v2[1000]; | |
584 | sprintf(v2,"%s>>hpullb",var); | |
585 | comp.fTree->Draw(v2, cut); | |
586 | TH1F * res = MakeCumul(hpullb); | |
587 | delete hpullb; | |
588 | return res; | |
589 | } | |
590 | ||
591 | ||
592 | ||
593 | void MakeComp(){ | |
594 | // | |
595 | // Set tree and set aliases | |
596 | // | |
597 | TFile *f = new TFile("trdComp1.root"); | |
598 | TTree * tree= (TTree*)f->Get("Comparison"); | |
599 | comp.fTree =tree; | |
600 | } | |
601 | ||
602 | void MakeAlias(){ | |
603 | ||
604 | comp.fTree->SetAlias("Pt","sqrt(MC.fParticle.fPx**2+MC.fParticle.fPy**2)"); | |
605 | comp.fTree->SetAlias("P","sqrt(MC.fParticle.fPx**2+MC.fParticle.fPy**2+MC.fParticle.fPz**2)"); | |
606 | comp.fTree->SetAlias("Ptpc","sqrt(fTrackRef.fPx**2+fTrackRef.fPy**2+fTrackRef.fPz**2)"); | |
607 | comp.fTree->SetAlias("VRadius","sqrt(MC.fParticle.fVx**2+MC.fParticle.fVy**2)"); | |
608 | comp.fTree->SetAlias("DRadius","sqrt(MC.fTRdecay.fX**2+MC.fTRdecay.fY**2)"); | |
609 | comp.fTree->SetAlias("ThetaM","MC.fParticle.fPz/Pt"); | |
610 | comp.fTree->SetAlias("rel0","(MC.fTOFReferences[0].P()-MC.fTPCReferences[0].P())/MC.fTOFReferences[0].P()"); //relative en. loss | |
611 | comp.fTree->SetAlias("rel1","(TR.P()-MC.fTPCReferences[0].P())/MC.fTPCReferences[0].P()"); //relative en. loss - in the last TRD | |
612 | comp.fTree->SetAlias("pt_tof","MC.fTOFReferences[0].Pt()"); | |
613 | comp.fTree->SetAlias("dpt_tof","(abs(1/RC.fRp[4])-MC.fTOFReferences[0].Pt())/MC.fTOFReferences[0].Pt()"); | |
614 | ||
615 | ||
616 | comp.fTree->SetAlias("d1pt_tof","(abs(RC.fRp[4])-1/MC.fTOFReferences[0].Pt())"); | |
617 | comp.fTree->SetAlias("d1ptp_tof","(abs(RC.fRp[4])-1/MC.fTOFReferences[0].Pt())/sqrt(RC.fRc[14])"); | |
618 | ||
619 | comp.fTree->SetAlias("dpt_trd","(abs(1/RC.fRp[4])-TR.Pt())/TR.Pt()"); | |
620 | ||
621 | comp.fTree->SetAlias("tofth","(MC.fTOFReferences[0].fPz/MC.fTOFReferences[0].Pt())"); | |
622 | comp.fTree->SetAlias("dth_tof","TRD0.fT[3]-tofth"); | |
7dd16671 | 623 | comp.fTree->SetAlias("pullth0","(TRD0.fT[3]-tofth)/sqrt(TRD0.fCtt)"); |
9dc17dc7 | 624 | comp.fTree->SetAlias("phi_tof","atan2(MC.fTOFReferences[0].fPy,MC.fTOFReferences[0].fPx)"); |
625 | comp.fTree->SetAlias("phi_trd","asin(TRD0.GetSnp())+TRD0.fAlpha"); | |
7dd16671 | 626 | comp.fTree->SetAlias("pullphi0","(sin(asin(TRD0.GetSnp())+TRD0.fAlpha)-sin(atan2(MC.fTOFReferences[0].fPy,MC.fTOFReferences[0].fPx)))/sqrt(TRD0.fCee)"); |
627 | comp.fTree->SetAlias("dphir0","(asin(TRD0.GetSnp())+TRD0.fAlpha-atan2(MC.fTOFReferences[0].fPy,MC.fTOFReferences[0].fPx))"); | |
9dc17dc7 | 628 | comp.fTree->SetAlias("dz0","(TRD0.fZ-MC.fTOFReferences[0].fZ)"); |
629 | comp.fTree->SetAlias("pully0","TRD0.fY/sqrt(TRD0.fCyy)"); | |
630 | comp.fTree->SetAlias("pullz0","(TRD0.fZ-MC.fTOFReferences[0].fZ)/sqrt(TRD0.fCzz)"); | |
631 | // | |
632 | comp.fTree->SetAlias("sigmay","sqrt(TRD0.fCyy)*(1+1/(1+abs(4./(RC.fRp[4]))))"); | |
633 | comp.fTree->SetAlias("pully3","TRD0.fY/sigmay"); | |
634 | comp.fTree->SetAlias("deltaz0","(TRD0.fZ-MC.fTOFReferences[0].fZ)"); | |
635 | // | |
636 | comp.fTree->SetAlias("pullc1","(abs(RC.fRp[4])-1/TR.Pt())/sqrt(RC.fRc[14])"); | |
637 | comp.fTree->SetAlias("deltac1","(abs(RC.fRp[4])-1/TR.Pt())*TR.Pt()"); | |
638 | comp.fTree->SetAlias("dist","sqrt(deltaz0**2+TRD0.fY**2)"); | |
639 | comp.fTree->SetAlias("pull0","pully0**2+pullz0**2"); | |
640 | comp.fTree->SetAlias("pullc0","(abs(TRD0.fC)-(1/(TRD0.fLocalConvConst*MC.fTOFReferences[0].Pt())))/sqrt(TRD0.fCcc)"); | |
641 | // | |
642 | comp.fTree->SetAlias("prob5","exp(-(TRD1.fTracklets[5].fChi2)/16.)*max((TRD1.fTracklets[5].fNFound-6)/9.,0)"); | |
643 | comp.fTree->SetAlias("prob4","exp(-(TRD1.fTracklets[4].fChi2)/16.)*max((TRD1.fTracklets[4].fNFound-6)/9.,0)"); | |
644 | // | |
e2cd3644 | 645 | comp.fTree->SetAlias("TPCE","sqrt(MC.fTPCReferences[2].P()**2+MC.fMass**2)"); |
9dc17dc7 | 646 | comp.fTree->SetAlias("TOFE","sqrt(MC.fTOFReferences[0].P()**2+MC.fMass**2)"); |
647 | comp.fTree->SetAlias("TRDE","sqrt(TR.P()**2+MC.fMass**2)"); | |
648 | comp.fTree->SetAlias("dtpi","TRD0.fIntegratedTime[2]-10^12*MC.fTOFReferences.fTime"); | |
649 | comp.fTree->SetAlias("dtk","TRD0.fIntegratedTime[3]-10^12*MC.fTOFReferences.fTime"); | |
650 | comp.fTree->SetAlias("dtp","TRD0.fIntegratedTime[4]-10^12*MC.fTOFReferences.fTime"); | |
651 | comp.fTree->SetAlias("dt","(abs(MC.fPdg)==2212)*dtp+(abs(MC.fPdg)==211)*dtpi+(abs(MC.fPdg)==321)*dtk"); //delta time | |
652 | comp.fTree->SetAlias("Beta","sqrt(TR.P()**2/(TR.P()**2+MC.fMass**2))"); | |
653 | comp.fTree->SetAlias("dtsigma","10+(1-Beta)*280"); | |
654 | comp.fTree->SetAlias("dptrel","(abs(TRD0.GetPt())-MC.fTOFReferences[0].Pt())/MC.fTOFReferences[0].Pt()"); | |
655 | // | |
656 | comp.fTree->SetAlias("dphi","atan2(TR.fY,TR..fX)-atan2(TR.fPy,TR..fPx)"); | |
657 | comp.fTree->SetAlias("dphi0","acos((TR.fX*TR.fPx+TR.fY*TR.fPy)/(TR.Pt()*TR.R()))"); | |
7dd16671 | 658 | comp.fTree->SetAlias("derel","(TRD0.fDE-(TPCE-TOFE))/(TPCE-TOFE)"); |
659 | // | |
660 | // TOF stuff | |
661 | // | |
662 | comp.fTree->SetAlias("tofm0","InfoTOF.fCl0.fLab[0]==abs(MC.fLabel)"); | |
663 | comp.fTree->SetAlias("tofm1","InfoTOF.fCl1.fLab[0]==abs(MC.fLabel)"); | |
664 | comp.fTree->SetAlias("tofm0s","InfoTOF.fCl0.fLab[0]>=abs(MC.fParticle.fDaughter[0]) && InfoTOF.fCl0.fLab[0]<=abs(MC.fParticle.fDaughter[1])"); | |
665 | comp.fTree->SetAlias("tofm1s","InfoTOF.fCl1.fLab[0]>=abs(MC.fParticle.fDaughter[0]) && InfoTOF.fCl1.fLab[0]<=abs(MC.fParticle.fDaughter[1])"); | |
666 | comp.fTree->SetAlias("tofm","tofm0||tofm1||tofm0s||tofm1s"); | |
667 | ||
668 | comp.fTree->SetAlias("Yref3","InfoMC.fSeeds[3].fYref[0]+InfoMC.fSeeds[3].fYref[1]*(InfoMC.fRef[3].LocalX()-InfoMC.fSeeds[3].fX0)"); | |
669 | comp.fTree->SetAlias("dyref3","Yref3-InfoMC.fRef[3].LocalY()"); | |
9dc17dc7 | 670 | |
7dd16671 | 671 | comp.fTree->SetAlias("Yref5","InfoMC.fSeeds[5].fYref[0]+InfoMC.fSeeds[5].fYref[1]*(InfoMC.fRef[5].LocalX()-InfoMC.fSeeds[5].fX0)"); |
672 | comp.fTree->SetAlias("dyref5","Yref5-InfoMC.fRef[5].LocalY()"); | |
9dc17dc7 | 673 | |
674 | } | |
675 | ||
676 | ||
677 | void MakeCompTr(){ | |
678 | // | |
679 | // Set tree and set aliases | |
680 | // | |
681 | TFile *ff = new TFile("TRDdebug.root"); | |
682 | TTree * treetr= (TTree*)ff->Get("tracklet"); | |
683 | comptr.fTree =treetr; | |
684 | ||
685 | } | |
686 | ||
687 | ||
688 | void DrawDE(TCut cut){ | |
689 | comp.DrawXY("(TRD0.fDE)","TPCE-TOFE",c2,cut,5,0.01,0.1,0,0.2); | |
690 | comp.fRes->Draw(); | |
691 | } | |
692 | ||
693 | ||
7dd16671 | 694 | |
695 | ||
9dc17dc7 | 696 | void DrawYResol1(TCut cut){ |
697 | comp.fTree->Draw("TRD1.fY","MC.fNTOFRef>0&&MC.fNTRDRef>5"+c1+cut,""); | |
698 | } | |
699 | ||
700 | void DrawYResol0(TCut cut){ | |
701 | comp.fTree->Draw("TRD0.fY","MC.fNTOFRef>0&&MC.fNTRDRef>5"+c0+cut,""); | |
702 | } | |
703 | ||
704 | ||
e2cd3644 | 705 | void DrawTOFResY(Int_t ndiv =10,Float_t ptmin=0.5, Float_t ptmax =1.5, Float_t dymin=-0.5, Float_t dymax=0.5){ |
9dc17dc7 | 706 | // |
7dd16671 | 707 | // |
708 | //Float_t ptmin=0.5, ptmax =1.5, dymin=-0.5, dymax=0.5; | |
709 | //Float_t ptmin=0.5, ptmax =5, dymin=-0.3, dymax=0.3; | |
e2cd3644 | 710 | // Int_t ndiv =5; |
7dd16671 | 711 | TCut cl("cl","InfoCl.fNClusters>60"); |
712 | comp.DrawXY("TR.Pt()","TRD0.fY",c0+c2+cl,"1",ndiv,ptmin, ptmax,dymin, dymax,60); | |
9dc17dc7 | 713 | TH1F * histof = (TH1F*)comp.fRes->Clone(); |
7dd16671 | 714 | comp.DrawXY("TR.Pt()","TRD1.fY",c0+c2+cl,"1",ndiv,ptmin, ptmax,dymin, dymax,60); |
9dc17dc7 | 715 | TH1F * histrd = (TH1F*)comp.fRes->Clone(); |
716 | histof->SetMarkerStyle(24); | |
717 | histrd->SetMarkerStyle(25); | |
718 | histof->SetXTitle("P_{t}(GeV)"); | |
719 | histof->SetYTitle("#Delta_{r-#phi}(cm)"); | |
720 | histof->Draw(); | |
721 | histrd->Draw("same"); | |
722 | TLegend *legend = new TLegend(0.55,0.12,0.85,0.35, "TRD resolution r-#phi (cm)"); | |
723 | legend->SetBorderSize(1); | |
724 | legend->AddEntry(histof ,"Resolution in TOF plane"); | |
725 | legend->AddEntry(histrd ,"Resolution in last TRD plane"); | |
726 | legend->Draw(); | |
727 | } | |
728 | ||
7dd16671 | 729 | void DrawTOFResZ(Float_t ptmin=0.5, Float_t ptmax =1.5, Float_t dzmin=-2.5, Float_t dzmax=2.5){ |
9dc17dc7 | 730 | // |
7dd16671 | 731 | comp.DrawXY("TR.Pt()","TRD0.fZ-MC.fTOFReferences.fZ",c0+c2,"1",10,ptmin,ptmax,dzmin,dzmax,30); |
9dc17dc7 | 732 | TH1F * histof = (TH1F*)comp.fRes->Clone(); |
7dd16671 | 733 | comp.DrawXY("TR.Pt()","TRD1.fZ-TR.fZ",c0+c2,"1",10,ptmin,ptmax,dzmin,dzmax,30); |
9dc17dc7 | 734 | TH1F * histrd = (TH1F*)comp.fRes->Clone(); |
735 | histof->SetMarkerStyle(24); | |
736 | histrd->SetMarkerStyle(25); | |
737 | histof->SetXTitle("P_{t}(GeV)"); | |
738 | histof->SetYTitle("#Delta_{z}(cm)"); | |
739 | histof->Draw(); | |
740 | histrd->Draw("same"); | |
741 | TLegend *legend = new TLegend(0.55,0.12,0.85,0.35, "TRD resolution z (cm)"); | |
742 | legend->SetBorderSize(1); | |
743 | legend->AddEntry(histof ,"Resolution in TOF plane"); | |
744 | legend->AddEntry(histrd ,"Resolution in last TRD plane"); | |
745 | legend->Draw(); | |
746 | } | |
747 | ||
748 | ||
749 | void MakeCumul(TCut & cut) | |
750 | { | |
7dd16671 | 751 | // TCut cut ="1&&TR.Pt()>0.4"; |
9dc17dc7 | 752 | TH1F * hisy = PullCumul("abs(pully0)",cut+cprim+c0+c2+c3,1000,20); |
753 | TH1F * hisz = PullCumul("abs(pullz0)",cut+cprim+c0+c2+c3,1000,20); | |
754 | TH1F * hisyz = PullCumul("sqrt(pullz0**2+pully0**2)",cut+cprim+c0+c2+c3,1000,20); | |
755 | hisz->SetXTitle("Pull (unit)"); | |
756 | hisz->SetYTitle("Cumulative function"); | |
7dd16671 | 757 | hisy->SetLineStyle(1);hisz->SetLineStyle(2);hisyz->SetLineStyle(3); |
758 | hisz->Draw(); hisy->Draw("same"); | |
759 | //hisyz->Draw("same"); | |
9dc17dc7 | 760 | TLegend *legend = new TLegend(0.55,0.12,0.85,0.35, "TRD pull (unit)"); |
761 | legend->SetBorderSize(1); | |
762 | legend->AddEntry(hisy ,"Pull y"); | |
763 | legend->AddEntry(hisz,"Pull z"); | |
7dd16671 | 764 | // legend->AddEntry(hisyz,"Pull yz"); |
9dc17dc7 | 765 | legend->Draw(); |
766 | } | |
767 | ||
768 | void MakeCumulAbs(TCut & cut) | |
769 | { | |
7dd16671 | 770 | // TCut cut = "1"; |
9dc17dc7 | 771 | TH1F * hisy = PullCumul("abs(TRD0.fY)",cut+cprim+c0+c2+c3,1000,6); |
772 | TH1F * hisz = PullCumul("abs(deltaz0)",cut+cprim+c0+c2+c3,1000,6); | |
773 | TH1F * hisyz = PullCumul("sqrt(TRD0.fY**2+deltaz0**2)",cut+cprim+c0+c2+c3,1000,6); | |
774 | hisz->SetXTitle("Delta (cm)"); | |
775 | hisz->SetYTitle("Cumulative function"); | |
776 | // | |
7dd16671 | 777 | hisy->SetLineStyle(1);hisz->SetLineStyle(2);hisyz->SetLineStyle(3); |
9dc17dc7 | 778 | hisz->Draw(); hisy->Draw("same");hisyz->Draw("same"); |
779 | TLegend *legend = new TLegend(0.55,0.12,0.85,0.35, "TRD residulas (cm)"); | |
780 | legend->SetBorderSize(1); | |
781 | legend->AddEntry(hisy ,"Delta y (cm)"); | |
782 | legend->AddEntry(hisz,"Delta z (cm)"); | |
783 | legend->AddEntry(hisyz,"Delta (cm)"); | |
784 | legend->Draw(); | |
785 | } | |
786 | ||
787 | void LoadClusters(TObjArray *arr, Int_t event) | |
788 | { | |
789 | // | |
790 | // | |
791 | // event=0 | |
792 | printf("Load cluster - event%d\n",event); | |
793 | fLoader->SetEventNumber(event); | |
794 | AliLoader * trdl = (AliLoader*)fLoader->GetLoader("TRDLoader"); | |
795 | trdl->LoadRecPoints(); | |
796 | ||
797 | TTree * ClusterTree = trdl->TreeR(); | |
798 | if (!ClusterTree) return; | |
799 | TBranch * branch = ClusterTree->GetBranch("TRDcluster"); | |
800 | TObjArray *clusterArray = new TObjArray(1000); | |
801 | branch->SetAddress(&clusterArray); | |
802 | Int_t nEntries = (Int_t) ClusterTree->GetEntries(); | |
803 | // | |
804 | for (Int_t iEntry = 0; iEntry < nEntries; iEntry++) { | |
805 | ClusterTree->GetEvent(iEntry); | |
806 | for (Int_t icl =0; icl<clusterArray->GetEntriesFast();icl++){ | |
807 | AliTRDcluster *cl = (AliTRDcluster*)clusterArray->UncheckedAt(icl); | |
808 | if (!cl) continue; | |
9dc17dc7 | 809 | Int_t detector=cl->GetDetector(); |
810 | Int_t localTimeBin=cl->GetLocalTimeBin(); | |
811 | // Int_t sector=geom.GetSector(detector); | |
812 | Int_t plane=geom.GetPlane(detector); | |
7dd16671 | 813 | for (Int_t ilab=0; ilab<3; ilab++){ |
814 | Int_t label = cl->GetLabel(ilab); | |
815 | Int_t pos = label; | |
816 | if (label<0) continue; | |
817 | if (localTimeBin>21||localTimeBin<0) continue; | |
818 | if (plane<0||plane>5) continue; | |
819 | ||
820 | AliTRDInfoCl * info =( AliTRDInfoCl*)arr->At(pos); | |
821 | if (!info) { | |
822 | info = new AliTRDInfoCl; | |
823 | arr->AddAt(info,pos); | |
824 | } | |
825 | // | |
826 | if (info->fCount[plane][localTimeBin]>0){ | |
827 | if (TMath::Abs(info->fCl[plane][localTimeBin].GetY()-cl->GetY())<20 && | |
828 | TMath::Abs(info->fCl[plane][localTimeBin].GetZ()-cl->GetZ())<20){ | |
829 | info->fCount[plane][localTimeBin]++; | |
830 | } | |
831 | }else{ | |
832 | info->fCl[plane][localTimeBin]=*cl; | |
833 | info->fCount[plane][localTimeBin]++; | |
834 | } | |
835 | } | |
9dc17dc7 | 836 | } |
837 | } | |
7dd16671 | 838 | for (Int_t i=0;i<arr->GetEntriesFast();i++){ |
839 | AliTRDInfoCl * info =( AliTRDInfoCl*)arr->At(i); | |
840 | if (info) info->Update(); | |
841 | } | |
9dc17dc7 | 842 | } |
843 | ||
844 | ||
7dd16671 | 845 | |
846 | void ReadSeeds(Int_t eventNr, TObjArray *sarray) | |
847 | { | |
848 | // | |
849 | // read seeds form the file | |
850 | // | |
851 | AliTRDseed seeds[6]; | |
852 | AliTRDseed *pseeds[6]={&seeds[0],&seeds[1],&seeds[2],&seeds[3],&seeds[4],&seeds[5]}; | |
853 | Int_t label, label1, label2, seventNr, nused; | |
854 | Float_t fakeRatio; | |
855 | TFile f("TRDdebug.root"); | |
856 | TTree *tree =(TTree*)f.Get("Seeds2"); | |
857 | if (!tree) return; | |
858 | TBranch * br0 = tree->GetBranch("S0."); | |
859 | TBranch * br1 = tree->GetBranch("S1."); | |
860 | TBranch * br2 = tree->GetBranch("S2."); | |
861 | TBranch * br3 = tree->GetBranch("S3."); | |
862 | TBranch * br4 = tree->GetBranch("S4."); | |
863 | TBranch * br5 = tree->GetBranch("S5."); | |
864 | TBranch * brL = tree->GetBranch("Label"); | |
865 | TBranch * brL1 = tree->GetBranch("Label1"); | |
866 | TBranch * brL2 = tree->GetBranch("Label2"); | |
867 | TBranch * brE = tree->GetBranch("EventNr"); | |
868 | TBranch * brU = tree->GetBranch("NUsed"); | |
869 | TBranch * brFR = tree->GetBranch("FakeRatio"); | |
870 | brE->SetAddress(&seventNr); | |
871 | brL->SetAddress(&label); | |
872 | brL1->SetAddress(&label1); | |
873 | brL2->SetAddress(&label2); | |
874 | brU->SetAddress(&nused); | |
875 | brFR->SetAddress(&fakeRatio); | |
876 | br0->SetAddress(&pseeds[0]); | |
877 | br1->SetAddress(&pseeds[1]); | |
878 | br2->SetAddress(&pseeds[2]); | |
879 | br3->SetAddress(&pseeds[3]); | |
880 | br4->SetAddress(&pseeds[4]); | |
881 | br5->SetAddress(&pseeds[5]); | |
882 | // | |
883 | // delete all seeds | |
884 | for (Int_t ip=0;ip<sarray->GetEntriesFast();ip++){ | |
885 | AliTRDseed * seed = (AliTRDseed*)sarray->UncheckedAt(ip); | |
886 | if (seed){ | |
887 | delete seed; | |
888 | } | |
889 | } | |
890 | sarray->Clear(); | |
891 | for (Int_t ip=0;ip<tree->GetEntries();ip++){ | |
892 | tree->GetEntry(ip); | |
893 | if (seventNr!=eventNr) continue; | |
894 | //if (nused>15) continue; | |
895 | if (!sarray->At(label)){ | |
896 | AliTRDInfoSeed * seed = new AliTRDInfoSeed; | |
897 | seed->Update(seeds); | |
898 | seed->fLabel = label; | |
899 | seed->fLabel1 = label1; | |
900 | seed->fLabel2 = label2; | |
901 | seed->fFakeRatio = fakeRatio; | |
902 | seed->fNUsed = nused; | |
903 | sarray->AddAt(seed,label); | |
904 | } | |
905 | } | |
906 | } | |
907 | ||
908 | void ReadTOFs(Int_t eventNr, TObjArray *tofs) | |
909 | { | |
910 | // | |
911 | // tofs->Clear(); | |
912 | // AliTOFtrackInfo *infotof = new AliTOFtrackInfo; | |
913 | // TFile ftof("TOFdebug.root"); | |
914 | // TTree * treetof =(TTree*)ftof.Get("Info"); | |
915 | // if (!treetof) return; | |
916 | // TBranch * br30 = treetof->GetBranch("Info."); | |
917 | // br30->SetAddress(&infotof); | |
918 | // // | |
919 | // // read tof info | |
920 | // for (Int_t ip=0;ip<treetof->GetEntries();ip++){ | |
921 | // treetof->GetEntry(ip); | |
922 | // if (eventNr!=infotof->fEventNr) continue; | |
923 | // // if (eventNr>info->fEventNr) continue; | |
924 | // // | |
925 | // Int_t label = TMath::Abs(infotof->fLab); | |
926 | // if (label==0) continue; | |
927 | // // if (label>=kmaxlabel) continue; | |
928 | // if (!(tofs->At(label))){ | |
929 | // tofs->AddAt(new AliTOFtrackInfo(*infotof),label); | |
930 | // } | |
931 | // else{ | |
932 | // if (infotof->fLength0>100) | |
933 | // tofs->AddAt(new AliTOFtrackInfo(*infotof),label); | |
934 | // } | |
935 | // } | |
936 | } | |
937 | ||
938 | void ReadTracks(Int_t eventNr, TObjArray *esds, TObjArray *trds) | |
939 | { | |
940 | // | |
941 | // read debug tracks | |
942 | // | |
943 | AliESDtrack * esdp = new AliESDtrack; | |
944 | AliTRDtrack * trdp = new AliTRDtrack; | |
945 | Int_t teventNr; | |
946 | TFile f2("TRDdebug.root"); | |
947 | TTree * tree2 =(TTree*)f2.Get("Tracks"); | |
948 | if (!tree2) return; | |
949 | TBranch * br20 = tree2->GetBranch("ESD."); | |
950 | TBranch * br21 = tree2->GetBranch("trd."); | |
951 | TBranch * br22 = tree2->GetBranch("EventNr"); | |
952 | br20->SetAddress(&esdp); | |
953 | br21->SetAddress(&trdp); | |
954 | br22->SetAddress(&teventNr); | |
955 | // | |
956 | // delete all event | |
957 | for (Int_t ip=0;ip<esds->GetEntriesFast();ip++){ | |
958 | if (esds->At(ip)){ | |
959 | delete esds->At(ip); | |
960 | delete trds->At(ip); | |
961 | } | |
962 | } | |
963 | esds->Clear(); | |
964 | trds->Clear(); | |
965 | for (Int_t ip=0;ip<tree2->GetEntries();ip++){ | |
966 | tree2->GetEntry(ip); | |
967 | if (eventNr!=teventNr) continue; | |
968 | // | |
e2cd3644 | 969 | Int_t label = TMath::Abs(esdp->GetLabel()); |
7dd16671 | 970 | if (label==0) continue; |
971 | if (!(esds->At(label))){ | |
972 | esds->AddAt(new AliESDtrack(*esdp),label); | |
973 | trds->AddAt(new AliTRDtrack(*trdp),label); | |
974 | } else{ | |
975 | AliTRDtrack * track =(AliTRDtrack*)trds->At(label); | |
976 | if (track->GetNumberOfClusters() <trdp->GetNumberOfClusters()){ | |
977 | esds->AddAt(new AliESDtrack(*esdp),label); | |
978 | trds->AddAt(new AliTRDtrack(*trdp),label); | |
979 | } | |
980 | } | |
981 | } | |
982 | } | |
983 | ||
984 | ||
985 | ||
9dc17dc7 | 986 | void MakeTree() |
987 | { | |
988 | fLoader = AliRunLoader::Open("galice.root"); | |
989 | AliESDtrack * esd = new AliESDtrack; | |
7dd16671 | 990 | AliTRDtrack * ptrd = new AliTRDtrack; |
991 | AliTRDtrack * ptrd0 = new AliTRDtrack; | |
992 | AliTRDtrack * ptrd1 = new AliTRDtrack; | |
993 | AliTRDtrack * ptrdb = new AliTRDtrack; | |
994 | // | |
995 | AliTRDtrack * trd = ptrd; | |
996 | AliTRDtrack * trd0 = ptrd0; | |
997 | AliTRDtrack * trd1 = ptrd1; | |
998 | AliTRDtrack * trdb = ptrdb; | |
999 | // | |
9dc17dc7 | 1000 | AliMCInfo* info = new AliMCInfo; |
1001 | AliTRDInfoMC *trdinfo = new AliTRDInfoMC; | |
7dd16671 | 1002 | AliTRDInfoCl *infocl = new AliTRDInfoCl; |
1003 | AliTRDInfoSeed *infoSeed = new AliTRDInfoSeed; | |
1004 | // AliTOFtrackInfo *infotof = new AliTOFtrackInfo; | |
1005 | // AliTOFtrackInfo dummytof; | |
1006 | // dummytof.fLab=-1; | |
1007 | AliTRDInfoSeed dummySeed; | |
1008 | AliTRDtrack dummytrd; | |
1009 | AliTRDtrack dummytrd0; | |
1010 | AliTRDtrack dummytrd1; | |
1011 | AliTRDtrack dummytrdb; | |
1012 | AliESDtrack dummyesd; | |
1013 | dummySeed.fLabel =-1; | |
9dc17dc7 | 1014 | AliTrackReference *ref = new AliTrackReference; |
7dd16671 | 1015 | TObjArray * clarray = new TObjArray(6*1000000); |
1016 | const Int_t kmaxlabel=1500000; | |
9dc17dc7 | 1017 | // |
1018 | // MC part | |
1019 | // | |
1020 | TFile f("cmpESDTracks.root"); | |
1021 | TTree * tree1 = (TTree*) f.Get("ESDcmpTracks"); | |
1022 | TBranch * brmc = tree1->GetBranch("MC"); | |
1023 | brmc->SetAddress(&info); | |
1024 | // | |
1025 | // | |
7dd16671 | 1026 | TObjArray *esds = new TObjArray(kmaxlabel); |
1027 | TObjArray *trds = new TObjArray(kmaxlabel); | |
1028 | TObjArray *tofs = new TObjArray(kmaxlabel); | |
1029 | TObjArray *seeds = new TObjArray(kmaxlabel); | |
1030 | // | |
9dc17dc7 | 1031 | // |
9dc17dc7 | 1032 | // |
9dc17dc7 | 1033 | // |
1034 | // comparison part | |
1035 | // | |
1036 | TFile f3("trdComp1.root","new"); | |
1037 | TTree * tree3 = new TTree("Comparison","Comparison"); | |
1038 | tree3->Branch("MC.","AliMCInfo",&info); | |
1039 | tree3->Branch("RC.","AliESDtrack",&esd); | |
1040 | tree3->Branch("TRD.","AliTRDtrack",&trd); | |
1041 | tree3->Branch("TRD0.","AliTRDtrack",&trd0); // track at tof ref | |
1042 | tree3->Branch("TRD1.","AliTRDtrack",&trd1); // track at last TRD ref | |
1043 | tree3->Branch("TRDb.","AliTRDtrack",&trdb); // track at tof ref | |
1044 | tree3->Branch("TR.","AliTrackReference",&ref); // track at last TRD ref | |
1045 | tree3->Branch("InfoMC.","AliTRDInfoMC",&trdinfo); // info about TRD track MC | |
1046 | tree3->Branch("InfoCl.","AliTRDInfoCl",&infocl); // info about TRD track MC | |
7dd16671 | 1047 | tree3->Branch("InfoSeed.","AliTRDInfoSeed",&infoSeed); // info about TRD track MC |
1048 | // tree3->Branch("InfoTOF.","AliTOFtrackInfo",&infotof); // info about TOF | |
9dc17dc7 | 1049 | // |
1050 | // | |
1051 | Int_t lastevent=-1; | |
1052 | ||
1053 | for (Int_t i=0;i<tree1->GetEntries();i++){ | |
1054 | brmc->GetEntry(i); | |
7dd16671 | 1055 | // if (info->fEventNr>5) break; |
9dc17dc7 | 1056 | if (lastevent!=info->fEventNr){ |
1057 | printf("Read Event\t%d\n", info->fEventNr); | |
7dd16671 | 1058 | ReadSeeds(info->fEventNr, seeds); |
1059 | ReadTOFs(info->fEventNr, tofs); | |
1060 | ReadTracks(info->fEventNr,esds,trds); | |
9dc17dc7 | 1061 | lastevent = info->fEventNr; |
9dc17dc7 | 1062 | clarray->Clear(); |
1063 | LoadClusters(clarray,lastevent); | |
1064 | // | |
9dc17dc7 | 1065 | } |
9dc17dc7 | 1066 | // if (info->fNTOFRef==0 && info->fNTRDRef<5) continue; |
1067 | Int_t label = info->fLabel; | |
1068 | trdinfo->Update(info); | |
1069 | if (clarray->At(label)){ | |
1070 | infocl = (AliTRDInfoCl*)clarray->At(label); | |
7dd16671 | 1071 | infocl->Update(); |
1072 | } | |
1073 | if (seeds->At(label)){ | |
1074 | infoSeed = (AliTRDInfoSeed*)seeds->At(label); | |
1075 | }else{ | |
1076 | infoSeed = &dummySeed; | |
9dc17dc7 | 1077 | } |
7dd16671 | 1078 | // if (info->fNTPCRef<3) continue; |
9dc17dc7 | 1079 | trd0->SetStop(kTRUE); |
1080 | trdb->SetStop(kTRUE); | |
1081 | trd1->SetStop(kTRUE); | |
7dd16671 | 1082 | // infotof = (AliTOFtrackInfo*)(tofs->At(label)); |
1083 | // if (!infotof) { | |
1084 | // infotof = &dummytof; | |
1085 | // infotof->fNCluster=0; | |
1086 | // esd = (AliESDtrack*)(esds->At(label)); | |
1087 | // if (esd){ | |
1088 | // Double_t pid[5]; | |
1089 | // esd->GetTPCpid(pid); | |
1090 | // Double_t suml=0; | |
1091 | // for (Int_t ip=0;ip<5;ip++) suml+=pid[ip]; | |
1092 | // for (Int_t ip=0;ip<5;ip++) { | |
1093 | // if (suml>0.0000000001) infotof->fRefPID[ip]= pid[ip]/suml; | |
1094 | // else infotof->fRefPID[ip]=0.2; | |
1095 | // } | |
1096 | // } | |
1097 | // } | |
1098 | trd = &dummytrd; | |
1099 | trd0 = &dummytrd0; | |
1100 | trd1 = &dummytrd1; | |
1101 | trdb = &dummytrdb; | |
1102 | esd = &dummyesd; | |
9dc17dc7 | 1103 | // |
1104 | if (esds->At(label)) { | |
7dd16671 | 1105 | trd0 = ptrd0; |
1106 | trd1 = ptrd1; | |
1107 | trdb = ptrdb; | |
9dc17dc7 | 1108 | // if (i%1000==0) printf("\n%d\n",i); |
1109 | trd = (AliTRDtrack*)(trds->At(label)); | |
1110 | esd = (AliESDtrack*)(esds->At(label)); | |
1111 | // | |
1112 | *trd0 = *trd; | |
1113 | if (trd->GetBackupTrack()) | |
1114 | *trdb = *trd->GetBackupTrack(); | |
1115 | else | |
1116 | *trdb = *trd; | |
1117 | *trd1 = *trd; | |
1118 | if (info->fNTOFRef>0){ | |
1119 | AliTrackReference * ref = ((AliTrackReference*)(info->fTOFReferences->At(0))); | |
1120 | Double_t x = ((AliTrackReference*)(info->fTOFReferences->At(0)))->X(); | |
1121 | Double_t y = ((AliTrackReference*)(info->fTOFReferences->At(0)))->Y(); | |
1122 | Double_t phi = TMath::ATan2(y,x); | |
1123 | Double_t radius = TMath::Sqrt(x*x+y*y); | |
1124 | // | |
1125 | // rotate and propagate to local system | |
1126 | // | |
1127 | trd0->SetStop(kFALSE); | |
1128 | if (!trd0->Rotate(phi-trd0->GetAlpha())) trd0->SetStop(kTRUE); | |
1129 | else{ | |
1130 | Double_t xyz0[3], param[7]; | |
1131 | trd0->GetGlobalXYZ(xyz0[0],xyz0[1],xyz0[2]); | |
1132 | Double_t xyz1[3]={ref->X(),ref->Y(),ref->Z()}; | |
1133 | AliKalmanTrack::MeanMaterialBudget(xyz0,xyz1,param); | |
1134 | if (!trd0->PropagateTo(radius,param[1],param[0])) trd0->SetStop(kTRUE); | |
1135 | } | |
1136 | trdb->SetStop(kFALSE); | |
1137 | if (!trdb->Rotate(phi-trdb->GetAlpha())) trdb->SetStop(kTRUE); | |
1138 | else | |
1139 | if (!trdb->PropagateTo(radius)) trdb->SetStop(kTRUE); | |
9dc17dc7 | 1140 | } |
1141 | ||
1142 | if (info->fNTRDRef>0){ | |
1143 | Double_t cradius=0; | |
1144 | Double_t cphi=0; | |
1145 | for (Int_t itrd =0; itrd<info->fNTRDRef;itrd++){ | |
1146 | Double_t x = ((AliTrackReference*)(info->fTRDReferences->At(itrd)))->X(); | |
1147 | Double_t y = ((AliTrackReference*)(info->fTRDReferences->At(itrd)))->Y(); | |
1148 | Double_t phi = TMath::ATan2(y,x); | |
1149 | // Double_t snp = TMath::ATan2(((AliTrackReference*)(info->fTRDReferences->At(itrd)))->Py(), | |
1150 | // ((AliTrackReference*)(info->fTRDReferences->At(itrd)))->Px()); | |
1151 | // if (TMath::Abs(snp)>0.9 &&cradius>10) break; | |
1152 | Double_t radius = TMath::Sqrt(x*x+y*y); | |
1153 | if (radius>cradius){ | |
1154 | cradius = radius; | |
1155 | cphi = phi; | |
1156 | ref = ((AliTrackReference*)(info->fTRDReferences->At(itrd))); | |
1157 | }else break; | |
1158 | } | |
1159 | // | |
1160 | // rotate and propagate to local system of the last track reference | |
1161 | // | |
1162 | trd1->SetStop(kFALSE); | |
1163 | if (!trd1->Rotate(cphi-trd1->GetAlpha())) trd1->SetStop(kTRUE); | |
1164 | else{ | |
1165 | if (!trd1->PropagateTo(cradius)) trd1->SetStop(kTRUE); | |
1166 | } | |
1167 | } | |
9dc17dc7 | 1168 | } |
7dd16671 | 1169 | tree3->Fill(); |
9dc17dc7 | 1170 | } |
1171 | f3.cd(); | |
1172 | tree3->Write(); | |
1173 | } | |
1174 | ||
1175 |