]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/qaRec/AliTRDresolution.cxx
Major update for the TRD tracking code
[u/mrichter/AliRoot.git] / TRD / qaRec / AliTRDresolution.cxx
CommitLineData
77203477 1/**************************************************************************
d2381af5 2* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3* *
4* Author: The ALICE Off-line Project. *
5* Contributors are mentioned in the code where appropriate. *
6* *
7* Permission to use, copy, modify and distribute this software and its *
8* documentation strictly for non-commercialf purposes is hereby granted *
9* without fee, provided that the above copyright notice appears in all *
10* copies and that both the copyright notice and this permission notice *
11* appear in the supporting documentation. The authors make no claims *
12* about the suitability of this software for any purpose. It is *
13* provided "as is" without express or implied warranty. *
14**************************************************************************/
77203477 15
e3cf3d02 16/* $Id: AliTRDresolution.cxx 27496 2008-07-22 08:35:45Z cblume $ */
77203477 17
18////////////////////////////////////////////////////////////////////////////
19// //
017bd6af 20// TRD tracking resolution //
21//
22// The class performs resolution and residual studies
23// of the TRD tracks for the following quantities :
24// - spatial position (y, [z])
25// - angular (phi) tracklet
26// - momentum at the track level
27//
28// The class has to be used for regular detector performance checks using the official macros:
29// - $ALICE_ROOT/TRD/qaRec/run.C
30// - $ALICE_ROOT/TRD/qaRec/makeResults.C
31//
32// For stand alone usage please refer to the following example:
33// {
34// gSystem->Load("libANALYSIS.so");
35// gSystem->Load("libTRDqaRec.so");
e3cf3d02 36// AliTRDresolution *res = new AliTRDresolution();
017bd6af 37// //res->SetMCdata();
38// //res->SetVerbose();
39// //res->SetVisual();
40// res->Load("TRD.TaskResolution.root");
41// if(!res->PostProcess()) return;
42// res->GetRefFigure(0);
43// }
44//
77203477 45// Authors: //
017bd6af 46// Alexandru Bercuci <A.Bercuci@gsi.de> //
77203477 47// Markus Fasel <M.Fasel@gsi.de> //
48// //
49////////////////////////////////////////////////////////////////////////////
50
aaf47b30 51#include <cstring>
52
124d488a 53#include <TROOT.h>
017bd6af 54#include <TSystem.h>
a6264da2 55#include <TPDGCode.h>
77203477 56#include <TObjArray.h>
7102d1b1 57#include <TH2.h>
58#include <TH1.h>
59#include <TF1.h>
017bd6af 60#include <TCanvas.h>
b2dc316d 61#include <TBox.h>
77203477 62#include <TProfile.h>
7102d1b1 63#include <TGraphErrors.h>
77203477 64#include <TMath.h>
aaf47b30 65#include "TTreeStream.h"
66#include "TGeoManager.h"
77203477 67
68#include "AliAnalysisManager.h"
77203477 69#include "AliTrackReference.h"
aaf47b30 70#include "AliTrackPointArray.h"
71#include "AliCDBManager.h"
72
9605ce80 73#include "AliTRDSimParam.h"
74#include "AliTRDgeometry.h"
75#include "AliTRDpadPlane.h"
aaf47b30 76#include "AliTRDcluster.h"
77#include "AliTRDseedV1.h"
78#include "AliTRDtrackV1.h"
79#include "AliTRDtrackerV1.h"
80#include "AliTRDReconstructor.h"
81#include "AliTRDrecoParam.h"
77203477 82
b2dc316d 83#include "AliTRDtrackInfo/AliTRDclusterInfo.h"
77203477 84#include "AliTRDtrackInfo/AliTRDtrackInfo.h"
e3cf3d02 85#include "AliTRDresolution.h"
77203477 86
e3cf3d02 87ClassImp(AliTRDresolution)
77203477 88
89//________________________________________________________
e3cf3d02 90AliTRDresolution::AliTRDresolution()
91 :AliTRDrecoTask("Resolution", "Spatial and Momentum Resolution")
017bd6af 92 ,fStatus(0)
aaf47b30 93 ,fReconstructor(0x0)
9605ce80 94 ,fGeo(0x0)
b718144c 95 ,fGraphS(0x0)
96 ,fGraphM(0x0)
6fc46cba 97 ,fCl(0x0)
98 ,fTrklt(0x0)
99 ,fMCcl(0x0)
100 ,fMCtrklt(0x0)
77203477 101{
aaf47b30 102 fReconstructor = new AliTRDReconstructor();
103 fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
9605ce80 104 fGeo = new AliTRDgeometry();
b2dc316d 105
de520d8f 106 InitFunctorList();
b2dc316d 107
6fc46cba 108 DefineOutput(1, TObjArray::Class()); // cluster2track
109 DefineOutput(2, TObjArray::Class()); // tracklet2track
110 DefineOutput(3, TObjArray::Class()); // cluster2mc
111 DefineOutput(4, TObjArray::Class()); // tracklet2mc
77203477 112}
113
ed383798 114//________________________________________________________
e3cf3d02 115AliTRDresolution::~AliTRDresolution()
ed383798 116{
b718144c 117 if(fGraphS){fGraphS->Delete(); delete fGraphS;}
118 if(fGraphM){fGraphM->Delete(); delete fGraphM;}
9605ce80 119 delete fGeo;
ed383798 120 delete fReconstructor;
2c0cf367 121 if(gGeoManager) delete gGeoManager;
6fc46cba 122 if(fCl){fCl->Delete(); delete fCl;}
123 if(fTrklt){fTrklt->Delete(); delete fTrklt;}
124 if(fMCcl){fMCcl->Delete(); delete fMCcl;}
125 if(fMCtrklt){fMCtrklt->Delete(); delete fMCtrklt;}
ed383798 126}
127
77203477 128
129//________________________________________________________
e3cf3d02 130void AliTRDresolution::CreateOutputObjects()
77203477 131{
132 // spatial resolution
77203477 133 OpenFile(0, "RECREATE");
39779ce6 134
3d86166d 135 fContainer = Histos();
b2dc316d 136
6fc46cba 137 fCl = new TObjArray();
138 fCl->SetOwner(kTRUE);
139 fTrklt = new TObjArray();
140 fTrklt->SetOwner(kTRUE);
141 fMCcl = new TObjArray();
142 fMCcl->SetOwner(kTRUE);
143 fMCtrklt = new TObjArray();
144 fMCtrklt->SetOwner(kTRUE);
77203477 145}
146
b2dc316d 147//________________________________________________________
e3cf3d02 148void AliTRDresolution::Exec(Option_t *opt)
b2dc316d 149{
6fc46cba 150 fCl->Delete();
151 fTrklt->Delete();
152 fMCcl->Delete();
153 fMCtrklt->Delete();
b2dc316d 154
155 AliTRDrecoTask::Exec(opt);
156
6fc46cba 157 PostData(1, fCl);
158 PostData(2, fTrklt);
159 PostData(3, fMCcl);
160 PostData(4, fMCtrklt);
b2dc316d 161}
aaf47b30 162
de520d8f 163//________________________________________________________
e3cf3d02 164TH1* AliTRDresolution::PlotCluster(const AliTRDtrackV1 *track)
de520d8f 165{
74b2e03d 166 if(track) fTrack = track;
167 if(!fTrack){
168 AliWarning("No Track defined.");
169 return 0x0;
de520d8f 170 }
171 TH1 *h = 0x0;
b1957d3c 172 if(!(h = ((TH2I*)fContainer->At(kCluster)))){
de520d8f 173 AliWarning("No output histogram defined.");
174 return 0x0;
175 }
176
b1957d3c 177 Double_t cov[3];
de520d8f 178 Float_t x0, y0, z0, dy, dydx, dzdx;
179 AliTRDseedV1 *fTracklet = 0x0;
180 for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
181 if(!(fTracklet = fTrack->GetTracklet(ily))) continue;
182 if(!fTracklet->IsOK()) continue;
de520d8f 183 x0 = fTracklet->GetX0();
184
185 // retrive the track angle with the chamber
0b8bcca4 186 y0 = fTracklet->GetYref(0);
187 z0 = fTracklet->GetZref(0);
188 dydx = fTracklet->GetYref(1);
189 dzdx = fTracklet->GetZref(1);
b1957d3c 190 fTracklet->GetCovRef(cov);
251a1ae6 191 Float_t tilt = fTracklet->GetTilt();
de520d8f 192 AliTRDcluster *c = 0x0;
193 fTracklet->ResetClusterIter(kFALSE);
194 while((c = fTracklet->PrevCluster())){
0b8bcca4 195 Float_t xc = c->GetX();
196 Float_t yc = c->GetY();
197 Float_t zc = c->GetZ();
198 Float_t dx = x0 - xc;
199 Float_t yt = y0 - dx*dydx;
200 Float_t zt = z0 - dx*dzdx;
b1957d3c 201 yc -= tilt*(zc-zt); // tilt correction
202 dy = yt - yc;
0b8bcca4 203
b1957d3c 204 h->Fill(dydx, dy/TMath::Sqrt(cov[0] + c->GetSigmaY2()));
de520d8f 205
206 if(fDebugLevel>=1){
de520d8f 207 // Get z-position with respect to anode wire
834ac2c9 208 //AliTRDSimParam *simParam = AliTRDSimParam::Instance();
0b8bcca4 209 Int_t istk = fGeo->GetStack(c->GetDetector());
de520d8f 210 AliTRDpadPlane *pp = fGeo->GetPadPlane(ily, istk);
211 Float_t row0 = pp->GetRow0();
58897a75 212 Float_t d = row0 - zt + pp->GetAnodeWireOffset();
de520d8f 213 d -= ((Int_t)(2 * d)) / 2.0;
214 if (d > 0.25) d = 0.5 - d;
b2dc316d 215
94f34623 216/* AliTRDclusterInfo *clInfo = new AliTRDclusterInfo;
6fc46cba 217 fCl->Add(clInfo);
0b8bcca4 218 clInfo->SetCluster(c);
219 clInfo->SetGlobalPosition(yt, zt, dydx, dzdx);
220 clInfo->SetResolution(dy);
221 clInfo->SetAnisochronity(d);
222 clInfo->SetDriftLength(dx);
223 (*fDebugStream) << "ClusterResiduals"
224 <<"clInfo.=" << clInfo
94f34623 225 << "\n";*/
de520d8f 226 }
227 }
228 }
229 return h;
230}
231
251a1ae6 232
233//________________________________________________________
e3cf3d02 234TH1* AliTRDresolution::PlotTracklet(const AliTRDtrackV1 *track)
251a1ae6 235{
236 if(track) fTrack = track;
237 if(!fTrack){
238 AliWarning("No Track defined.");
239 return 0x0;
240 }
241 TH1 *h = 0x0;
b1957d3c 242 if(!(h = ((TH2I*)fContainer->At(kTrackletY)))){
251a1ae6 243 AliWarning("No output histogram defined.");
244 return 0x0;
245 }
246
b1957d3c 247 Double_t cov[3], covR[3];
e3cf3d02 248 Float_t x, y0, dx, dy, dydx;
b1957d3c 249 AliTRDseedV1 *fTracklet = 0x0;
250 for(Int_t il=AliTRDgeometry::kNlayer; il--;){
251 if(!(fTracklet = fTrack->GetTracklet(il))) continue;
252 if(!fTracklet->IsOK()) continue;
253 y0 = fTracklet->GetYref(0);
254 dydx = fTracklet->GetYref(1);
e3cf3d02 255 x = fTracklet->GetX();
256 dx = fTracklet->GetX0() - x;
257 dy = y0-dx*dydx - fTracklet->GetY();
258 fTracklet->GetCovAt(x, cov);
b1957d3c 259 fTracklet->GetCovRef(covR);
260 h->Fill(dydx, dy/TMath::Sqrt(cov[0] + covR[0]));
251a1ae6 261 }
262 return h;
263}
264
265//________________________________________________________
e3cf3d02 266TH1* AliTRDresolution::PlotTrackletPhi(const AliTRDtrackV1 *track)
251a1ae6 267{
268 if(track) fTrack = track;
269 if(!fTrack){
270 AliWarning("No Track defined.");
271 return 0x0;
272 }
273 TH1 *h = 0x0;
b1957d3c 274 if(!(h = ((TH2I*)fContainer->At(kTrackletPhi)))){
251a1ae6 275 AliWarning("No output histogram defined.");
276 return 0x0;
277 }
a37c3c70 278
bcdbe5e5 279 AliTRDseedV1 *tracklet = 0x0;
a37c3c70 280 for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++){
a37c3c70 281 if(!(tracklet = fTrack->GetTracklet(il))) continue;
282 if(!tracklet->IsOK()) continue;
bcdbe5e5 283 h->Fill(tracklet->GetYref(1), TMath::ATan(tracklet->GetYref(1))-TMath::ATan(tracklet->GetYfit(1)));
251a1ae6 284 }
285 return h;
286}
287
288
de520d8f 289//________________________________________________________
e3cf3d02 290TH1* AliTRDresolution::PlotResolution(const AliTRDtrackV1 *track)
de520d8f 291{
c0811145 292 if(!HasMCdata()){
74b2e03d 293 AliWarning("No MC defined. Results will not be available.");
de520d8f 294 return 0x0;
295 }
74b2e03d 296 if(track) fTrack = track;
297 if(!fTrack){
298 AliWarning("No Track defined.");
299 return 0x0;
de520d8f 300 }
301 TH1 *h = 0x0;
612ae7ed 302 UChar_t s;
de520d8f 303 Int_t pdg = fMC->GetPDG(), det=-1;
0b8bcca4 304 Int_t label = fMC->GetLabel();
e3cf3d02 305 Double_t x, y, z;
306 Float_t p, pt, x0, y0, z0, dx, dy, dz, dydx, dzdx;
307 Double_t covR[3];
b1957d3c 308
0e08101e 309 if(fDebugLevel>=1){
310 Double_t DX[12], DY[12], DZ[12], DPt[12], COV[12][15];
311 fMC->PropagateKalman(DX, DY, DZ, DPt, COV);
312 (*fDebugStream) << "MCkalman"
313 << "pdg=" << pdg
314 << "dx0=" << DX[0]
315 << "dx1=" << DX[1]
316 << "dx2=" << DX[2]
317 << "dy0=" << DY[0]
318 << "dy1=" << DY[1]
319 << "dy2=" << DY[2]
320 << "dz0=" << DZ[0]
321 << "dz1=" << DZ[1]
322 << "dz2=" << DZ[2]
323 << "dpt0=" << DPt[0]
324 << "dpt1=" << DPt[1]
325 << "dpt2=" << DPt[2]
326 << "\n";
327 }
328
de520d8f 329 AliTRDseedV1 *fTracklet = 0x0;
330 for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
b1957d3c 331 if(!(fTracklet = fTrack->GetTracklet(ily)) ||
332 !fTracklet->IsOK()) continue;
333
de520d8f 334 det = fTracklet->GetDetector();
335 x0 = fTracklet->GetX0();
b1957d3c 336 //radial shift with respect to the MC reference (radial position of the pad plane)
e3cf3d02 337 x= fTracklet->GetX();
338 dx = x0 - x;
fd40f855 339 if(!fMC->GetDirections(x0, y0, z0, dydx, dzdx, pt, s)) continue;
e3cf3d02 340 if(fDebugLevel>=1){
341 (*fDebugStream) << "MC"
342 << "det=" << det
343 << "pdg=" << pdg
344 << "pt=" << pt
345 << "x0=" << x0
346 << "y0=" << y0
347 << "z0=" << z0
348 << "dydx=" << dydx
349 << "dzdx=" << dzdx
350 << "\n";
351 }
b1957d3c 352 // MC track position at reference radial position
e3cf3d02 353 Float_t yt = y0 - dx*dydx;
354 Float_t zt = z0 - dx*dzdx;
6090b78a 355 p = pt*(1.+dzdx*dzdx); // pt -> p
b1957d3c 356
6090b78a 357 // add Kalman residuals for y, z and pt
e3cf3d02 358 Float_t yr = fTracklet->GetYref(0) - dx*fTracklet->GetYref(1);
b1957d3c 359 dy = yt - yr;
fd40f855 360 Float_t zr = fTracklet->GetZref(0) - dx*fTracklet->GetZref(1);
b1957d3c 361 dz = zt - zr;
6090b78a 362 Float_t tgl = fTracklet->GetTgl();
e3cf3d02 363 Float_t dpt = pt - fTracklet->GetMomentum()/(1.+tgl*tgl);
364 fTracklet->GetCovRef(covR);
6090b78a 365
b1957d3c 366 ((TH2I*)fContainer->At(kMCtrackY))->Fill(dydx, dy);
e3cf3d02 367 ((TH2I*)fContainer->At(kMCtrackYPull))->Fill(dydx, dy/TMath::Sqrt(covR[0]));
368 if(ily==0){
369 ((TH2I*)fContainer->At(kMCtrackZIn))->Fill(dzdx, dz);
370 ((TH2I*)fContainer->At(kMCtrackZInPull))->Fill(dzdx, dz/TMath::Sqrt(covR[2]));
371 } else if(ily==AliTRDgeometry::kNlayer-1) {
372 ((TH2I*)fContainer->At(kMCtrackZOut))->Fill(dzdx, dz);
373 ((TH2I*)fContainer->At(kMCtrackZOutPull))->Fill(dzdx, dz/TMath::Sqrt(covR[2]));
374 }
375 if(pdg!=kElectron && pdg!=kPositron) ((TH2I*)fContainer->At(kMCtrackPt))->Fill(1./pt, dpt);
b1957d3c 376 // Fill Debug stream for Kalman track
377 if(fDebugLevel>=1){
b1957d3c 378 Float_t dydxr = fTracklet->GetYref(1);
b1957d3c 379 (*fDebugStream) << "MCtrack"
e3cf3d02 380 << "x=" << x
381 << "dpt=" << dpt
b1957d3c 382 << "dy=" << dy
383 << "dz=" << dz
384 << "dydxr=" << dydxr
6090b78a 385 << "dzdxr=" << tgl
e3cf3d02 386 << "s2y=" << covR[0]
387 << "s2z=" << covR[2]
b1957d3c 388 << "\n";
389 }
de520d8f 390
391 // recalculate tracklet based on the MC info
392 AliTRDseedV1 tt(*fTracklet);
393 tt.SetZref(0, z0);
394 tt.SetZref(1, dzdx);
e3cf3d02 395 tt.Fit(kTRUE);
396 x= tt.GetX(); // the true one
397 dx = x0 - x;
398 yt = y0 - dx*dydx;
399 Bool_t rc = tt.IsRowCross();
400
b1957d3c 401 // add tracklet residuals for y and dydx
e3cf3d02 402 Float_t yf = tt.GetY();
403 dy = yt - yf; dz = 100.;
b1957d3c 404 Float_t dphi = (tt.GetYfit(1) - dydx);
405 dphi /= 1.- tt.GetYfit(1)*dydx;
e3cf3d02 406 Double_t s2y = tt.GetS2Y(), s2z = tt.GetS2Z();
407 if(!rc){
408 ((TH2I*)fContainer->At(kMCtrackletY))->Fill(dydx, dy);
409 ((TH2I*)fContainer->At(kMCtrackletYPull))->Fill(dydx, dy/TMath::Sqrt(s2y));
410 ((TH2I*)fContainer->At(kMCtrackletPhi))->Fill(dydx, dphi*TMath::RadToDeg());
411 } else {
b1957d3c 412 // add tracklet residuals for z
e3cf3d02 413 dz = tt.GetZ() - (z0 - dx*dzdx) ;
b1957d3c 414 ((TH2I*)fContainer->At(kMCtrackletZ))->Fill(dzdx, dz);
e3cf3d02 415 ((TH2I*)fContainer->At(kMCtrackletZPull))->Fill(dzdx, dz/TMath::Sqrt(s2z));
de520d8f 416 }
417
b1957d3c 418 // Fill Debug stream for tracklet
de520d8f 419 if(fDebugLevel>=1){
b1957d3c 420 (*fDebugStream) << "MCtracklet"
e3cf3d02 421 << "rc=" << rc
422 << "x=" << x
423 << "dy=" << dy
424 << "dz=" << dz
425 << "dphi=" << dphi
426 << "s2y=" << s2y
427 << "s2z=" << s2z
de520d8f 428 << "\n";
429 }
39779ce6 430
de520d8f 431 Int_t istk = AliTRDgeometry::GetStack(det);
432 AliTRDpadPlane *pp = fGeo->GetPadPlane(ily, istk);
58897a75 433 Float_t zr0 = pp->GetRow0() + pp->GetAnodeWireOffset();
de520d8f 434 Float_t tilt = fTracklet->GetTilt();
435
436 AliTRDcluster *c = 0x0;
437 fTracklet->ResetClusterIter(kFALSE);
438 while((c = fTracklet->PrevCluster())){
439 Float_t q = TMath::Abs(c->GetQ());
b1957d3c 440 //AliTRDseedV1::GetClusterXY(c,x,y);
e3cf3d02 441 x = c->GetX(); y = c->GetY(); z = c->GetZ();
d937ad7a 442 Float_t xc = x;
443 Float_t yc = y;
e3cf3d02 444 Float_t zc = z;
a37c3c70 445 dx = x0 - xc;
fd40f855 446 yt = y0 - dx*dydx;
447 zt = z0 - dx*dzdx;
de520d8f 448 dy = yt - (yc - tilt*(zc-zt));
449
450 // Fill Histograms
b1957d3c 451 if(q>20. && q<250.) ((TH2I*)fContainer->At(kMCcluster))->Fill(dydx, dy);
452
453 // Fill calibration container
454 Float_t d = zr0 - zt;
455 d -= ((Int_t)(2 * d)) / 2.0;
456 if (d > 0.25) d = 0.5 - d;
457 AliTRDclusterInfo *clInfo = new AliTRDclusterInfo;
6fc46cba 458 fMCcl->Add(clInfo);
b1957d3c 459 clInfo->SetCluster(c);
460 clInfo->SetMC(pdg, label);
461 clInfo->SetGlobalPosition(yt, zt, dydx, dzdx);
462 clInfo->SetResolution(dy);
463 clInfo->SetAnisochronity(d);
fd40f855 464 clInfo->SetDriftLength(dx-.5*AliTRDgeometry::CamHght());
b1957d3c 465 clInfo->SetTilt(tilt);
b2dc316d 466
b1957d3c 467 // Fill Debug Tree
468 if(fDebugLevel>=2){
b2dc316d 469 //clInfo->Print();
b1957d3c 470 (*fDebugStream) << "MCcluster"
b2dc316d 471 <<"clInfo.=" << clInfo
de520d8f 472 << "\n";
473 }
474 }
77203477 475 }
de520d8f 476 return h;
77203477 477}
478
de520d8f 479
d85cd79c 480//________________________________________________________
e3cf3d02 481Bool_t AliTRDresolution::GetRefFigure(Int_t ifig)
d85cd79c 482{
b1957d3c 483 Float_t y[2] = {0., 0.};
b2dc316d 484 TBox *b = 0x0;
a391a274 485 TAxis *ax = 0x0;
486 TGraphErrors *g = 0x0;
d85cd79c 487 switch(ifig){
b1957d3c 488 case kCluster:
de520d8f 489 if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
a391a274 490 g->Draw("apl");
491 ax = g->GetHistogram()->GetYaxis();
b1957d3c 492 y[0] = -0.5; y[1] = 2.5;
493 ax->SetRangeUser(y[0], y[1]);
494 ax->SetTitle("Cluster-Track Pools #sigma/#mu [mm]");
a391a274 495 ax = g->GetHistogram()->GetXaxis();
017bd6af 496 ax->SetTitle("tg(#phi)");
de520d8f 497 if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
017bd6af 498 g->Draw("pl");
b1957d3c 499 b = new TBox(-.15, y[0], .15, y[1]);
b2dc316d 500 b->SetFillStyle(3002);b->SetFillColor(kGreen);
501 b->SetLineColor(0); b->Draw();
e15179be 502 return kTRUE;
b1957d3c 503 case kTrackletY:
251a1ae6 504 if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
505 g->Draw("apl");
506 ax = g->GetHistogram()->GetYaxis();
a37c3c70 507 ax->SetRangeUser(-.5, 3.);
b1957d3c 508 ax->SetTitle("Tracklet-Track Y-Pools #sigma/#mu [mm]");
251a1ae6 509 ax = g->GetHistogram()->GetXaxis();
510 ax->SetTitle("tg(#phi)");
511 if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
512 g->Draw("pl");
b1957d3c 513 b = new TBox(-.15, y[0], .15, y[1]);
251a1ae6 514 b->SetFillStyle(3002);b->SetFillColor(kGreen);
515 b->SetLineColor(0); b->Draw();
e15179be 516 return kTRUE;
b1957d3c 517 case kTrackletPhi:
251a1ae6 518 if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
519 g->Draw("apl");
520 ax = g->GetHistogram()->GetYaxis();
a37c3c70 521 ax->SetRangeUser(-.5, 2.);
b1957d3c 522 ax->SetTitle("Tracklet-Track Phi-Residuals #sigma/#mu [rad]");
251a1ae6 523 ax = g->GetHistogram()->GetXaxis();
524 ax->SetTitle("tg(#phi)");
525 if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
526 g->Draw("pl");
b1957d3c 527 b = new TBox(-.15, y[0], .15, y[1]);
251a1ae6 528 b->SetFillStyle(3002);b->SetFillColor(kGreen);
529 b->SetLineColor(0); b->Draw();
e15179be 530 return kTRUE;
b1957d3c 531 case kMCcluster:
de520d8f 532 if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
a391a274 533 ax = g->GetHistogram()->GetYaxis();
834ac2c9 534 y[0] = -.05; y[1] = 0.6;
b1957d3c 535 ax->SetRangeUser(y[0], y[1]);
536 ax->SetTitle("Y_{cluster} #sigma/#mu [mm]");
b718144c 537 ax = g->GetHistogram()->GetXaxis();
834ac2c9 538 ax->SetRangeUser(-.3, .3);
017bd6af 539 ax->SetTitle("tg(#phi)");
b718144c 540 g->Draw("apl");
de520d8f 541 if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
017bd6af 542 g->Draw("pl");
b1957d3c 543 b = new TBox(-.15, y[0], .15, y[1]);
544 b->SetFillStyle(3002);b->SetFillColor(kBlue);
b2dc316d 545 b->SetLineColor(0); b->Draw();
e15179be 546 return kTRUE;
b1957d3c 547 case kMCtrackletY:
de520d8f 548 if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
b718144c 549 ax = g->GetHistogram()->GetYaxis();
834ac2c9 550 y[0] = -.05; y[1] = 0.25;
b1957d3c 551 ax->SetRangeUser(y[0], y[1]);
552 ax->SetTitle("Y_{tracklet} #sigma/#mu [mm]");
a391a274 553 ax = g->GetHistogram()->GetXaxis();
834ac2c9 554 ax->SetRangeUser(-.2, .2);
612ae7ed 555 ax->SetTitle("tg(#phi)");
de520d8f 556 g->Draw("apl");
557 if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
558 g->Draw("pl");
b1957d3c 559 b = new TBox(-.15, y[0], .15, y[1]);
560 b->SetFillStyle(3002);b->SetFillColor(kBlue);
561 b->SetLineColor(0); b->Draw();
e15179be 562 return kTRUE;
b1957d3c 563 case kMCtrackletZ:
de520d8f 564 if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
565 ax = g->GetHistogram()->GetYaxis();
566 ax->SetRangeUser(-.5, 1.);
b1957d3c 567 ax->SetTitle("Z_{tracklet} #sigma/#mu [mm]");
de520d8f 568 ax = g->GetHistogram()->GetXaxis();
612ae7ed 569 ax->SetTitle("tg(#theta)");
a391a274 570 g->Draw("apl");
de520d8f 571 if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
017bd6af 572 g->Draw("pl");
e15179be 573 return kTRUE;
b1957d3c 574 case kMCtrackletPhi:
de520d8f 575 if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
a391a274 576 ax = g->GetHistogram()->GetYaxis();
b1957d3c 577 y[0] = -.05; y[1] = .2;
578 ax->SetRangeUser(y[0], y[1]);
579 ax->SetTitle("#Phi_{tracklet} #sigma/#mu [deg]");
a391a274 580 ax = g->GetHistogram()->GetXaxis();
de520d8f 581 ax->SetTitle("tg(#phi)");
a391a274 582 g->Draw("apl");
de520d8f 583 if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
017bd6af 584 g->Draw("pl");
e15179be 585 return kTRUE;
b1957d3c 586 case kMCtrackY:
587 if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
588 ax = g->GetHistogram()->GetYaxis();
834ac2c9 589 y[0] = -.05; y[1] = 0.25;
b1957d3c 590 ax->SetRangeUser(y[0], y[1]);
591 ax->SetTitle("Y_{track} #sigma/#mu [mm]");
592 ax = g->GetHistogram()->GetXaxis();
834ac2c9 593 ax->SetRangeUser(-.2, .2);
b1957d3c 594 ax->SetTitle("tg(#phi)");
595 g->Draw("apl");
596 if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
597 g->Draw("pl");
598 b = new TBox(-.15, y[0], .15, y[1]);
599 b->SetFillStyle(3002);b->SetFillColor(kBlue);
600 b->SetLineColor(0); b->Draw();
d937ad7a 601 return kTRUE;
e3cf3d02 602 case kMCtrackZIn:
b1957d3c 603 if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
604 ax = g->GetHistogram()->GetYaxis();
605 ax->SetRangeUser(-.5, 2.);
606 ax->SetTitle("Z_{track} #sigma/#mu [mm]");
607 ax = g->GetHistogram()->GetXaxis();
608 ax->SetTitle("tg(#theta)");
609 g->Draw("apl");
610 if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
611 g->Draw("pl");
d937ad7a 612 return kTRUE;
6090b78a 613 case kMCtrackPt:
614 if(!(g = (TGraphErrors*)fGraphS->At(ifig))) break;
615 ax = g->GetHistogram()->GetYaxis();
616 ax->SetRangeUser(-.5, 2.);
617 ax->SetTitle("#epsilon_{P_{t}}^{track} / #mu [%]");
618 ax = g->GetHistogram()->GetXaxis();
619 ax->SetTitle("1/p_{t}");
620 g->Draw("apl");
621 if(!(g = (TGraphErrors*)fGraphM->At(ifig))) break;
622 g->Draw("pl");
623 return kTRUE;
d85cd79c 624 }
017bd6af 625 AliInfo(Form("Reference plot [%d] missing result", ifig));
e15179be 626 return kFALSE;
d85cd79c 627}
628
39779ce6 629
77203477 630//________________________________________________________
e3cf3d02 631Bool_t AliTRDresolution::PostProcess()
7102d1b1 632{
d85cd79c 633 //fContainer = dynamic_cast<TObjArray*>(GetOutputData(0));
634 if (!fContainer) {
635 Printf("ERROR: list not available");
636 return kFALSE;
3d86166d 637 }
b718144c 638 fNRefFigures = fContainer->GetEntriesFast();
b1957d3c 639 TGraphErrors *gm = 0x0, *gs = 0x0;
b718144c 640 if(!fGraphS){
641 fGraphS = new TObjArray(fNRefFigures);
642 fGraphS->SetOwner();
b1957d3c 643 for(Int_t ig=0; ig<fNRefFigures; ig++){
644 gs = new TGraphErrors();
645 gs->SetLineColor(kRed);
646 gs->SetMarkerStyle(23);
647 gs->SetMarkerColor(kRed);
648 gs->SetNameTitle(Form("s_%d", ig), "");
649 fGraphS->AddAt(gs, ig);
650 }
b718144c 651 }
652 if(!fGraphM){
653 fGraphM = new TObjArray(fNRefFigures);
654 fGraphM->SetOwner();
b1957d3c 655 for(Int_t ig=0; ig<fNRefFigures; ig++){
656 gm = new TGraphErrors();
657 gm->SetLineColor(kBlue);
658 gm->SetMarkerStyle(7);
659 gm->SetMarkerColor(kBlue);
660 gm->SetNameTitle(Form("m_%d", ig), "");
661 fGraphM->AddAt(gm, ig);
662 }
b718144c 663 }
7102d1b1 664
d2381af5 665 TH2I *h2 = 0x0;
666 TH1D *h = 0x0;
b718144c 667
668 // define models
d2381af5 669 TF1 f("f1", "gaus", -.5, .5);
b718144c 670
017bd6af 671 TF1 fb("fb", "[0]*exp(-0.5*((x-[1])/[2])**2)+[3]", -.5, .5);
672
b718144c 673 TF1 fc("fc", "[0]*exp(-0.5*((x-[1])/[2])**2)+[3]*exp(-0.5*((x-[4])/[5])**2)", -.5, .5);
874acced 674
017bd6af 675 TCanvas *c = 0x0;
676 if(IsVisual()) c = new TCanvas("c", Form("%s Visual", GetName()), 500, 500);
677 char opt[5];
678 sprintf(opt, "%c%c", IsVerbose() ? ' ' : 'Q', IsVisual() ? ' ': 'N');
679
680
874acced 681 //PROCESS RESIDUAL DISTRIBUTIONS
682
683 // Clusters residuals
b1957d3c 684 h2 = (TH2I *)(fContainer->At(kCluster));
685 gm = (TGraphErrors*)fGraphM->At(kCluster);
686 gs = (TGraphErrors*)fGraphS->At(kCluster);
874acced 687 for(Int_t ibin = 1; ibin <= h2->GetNbinsX(); ibin++){
688 Double_t phi = h2->GetXaxis()->GetBinCenter(ibin);
874acced 689 h = h2->ProjectionY("py", ibin, ibin);
b2dc316d 690 if(h->GetEntries()<100) continue;
8d61815d 691 AdjustF1(h, &f);
017bd6af 692
693 if(IsVisual()){c->cd(); c->SetLogy();}
8d61815d 694 h->Fit(&f, opt, "", -0.5, 0.5);
017bd6af 695 if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
696
b2dc316d 697 Int_t ip = gm->GetN();
698 gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
699 gm->SetPointError(ip, 0., 10.*f.GetParError(1));
700 gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
701 gs->SetPointError(ip, 0., 10.*f.GetParError(2));
874acced 702 }
d2381af5 703
251a1ae6 704 // Tracklet y residuals
b1957d3c 705 h2 = (TH2I *)(fContainer->At(kTrackletY));
706 gm = (TGraphErrors*)fGraphM->At(kTrackletY);
707 gs = (TGraphErrors*)fGraphS->At(kTrackletY);
251a1ae6 708 for(Int_t ibin = 1; ibin <= h2->GetNbinsX(); ibin++){
709 Double_t phi = h2->GetXaxis()->GetBinCenter(ibin);
710 h = h2->ProjectionY("py", ibin, ibin);
711 if(h->GetEntries()<100) continue;
712 AdjustF1(h, &f);
713
714 if(IsVisual()){c->cd(); c->SetLogy();}
715 h->Fit(&f, opt, "", -0.5, 0.5);
716 if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
717
718 Int_t ip = gm->GetN();
719 gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
720 gm->SetPointError(ip, 0., 10.*f.GetParError(1));
721 gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
722 gs->SetPointError(ip, 0., 10.*f.GetParError(2));
723 }
724
725 // Tracklet phi residuals
b1957d3c 726 h2 = (TH2I *)(fContainer->At(kTrackletPhi));
727 gm = (TGraphErrors*)fGraphM->At(kTrackletPhi);
728 gs = (TGraphErrors*)fGraphS->At(kTrackletPhi);
251a1ae6 729 for(Int_t ibin = 1; ibin <= h2->GetNbinsX(); ibin++){
730 Double_t phi = h2->GetXaxis()->GetBinCenter(ibin);
731 h = h2->ProjectionY("py", ibin, ibin);
732 if(h->GetEntries()<100) continue;
733 AdjustF1(h, &f);
734
735 if(IsVisual()){c->cd(); c->SetLogy();}
736 h->Fit(&f, opt, "", -0.5, 0.5);
737 if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
738
739 Int_t ip = gm->GetN();
740 gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
741 gm->SetPointError(ip, 0., 10.*f.GetParError(1));
742 gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
743 gs->SetPointError(ip, 0., 10.*f.GetParError(2));
744 }
745
b1957d3c 746 if(!HasMCdata()){
747 if(c) delete c;
748 return kTRUE;
749 }
874acced 750
b1957d3c 751 //PROCESS MC RESIDUAL DISTRIBUTIONS
752
753 // cluster y resolution
754 h2 = (TH2I*)fContainer->At(kMCcluster);
755 gm = (TGraphErrors*)fGraphM->At(kMCcluster);
756 gs = (TGraphErrors*)fGraphS->At(kMCcluster);
757 for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
758 h = h2->ProjectionY("py", iphi, iphi);
759 if(h->GetEntries()<100) continue;
760 AdjustF1(h, &f);
761
762 if(IsVisual()){c->cd(); c->SetLogy();}
763 h->Fit(&f, opt, "", -0.5, 0.5);
764 if(IsVerbose()){
765 printf("phi[%d] mean[%e] sigma[%e]\n\n", iphi, 10.*f.GetParameter(1), 10.*f.GetParameter(2));
d2381af5 766 }
b1957d3c 767 if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
768
769 Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
770 Int_t ip = gm->GetN();
771 gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
772 gm->SetPointError(ip, 0., 10.*f.GetParError(1));
773 gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
774 gs->SetPointError(ip, 0., 10.*f.GetParError(2));
775 }
776
777 // tracklet y resolution
778 h2 = (TH2I*)fContainer->At(kMCtrackletY);
779 gm = (TGraphErrors*)fGraphM->At(kMCtrackletY);
780 gs = (TGraphErrors*)fGraphS->At(kMCtrackletY);
781 for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
782 h = h2->ProjectionY("py", iphi, iphi);
783 if(h->GetEntries()<100) continue;
784 AdjustF1(h, &f);
785
786 if(IsVisual()){c->cd(); c->SetLogy();}
787 h->Fit(&f, opt, "", -0.5, 0.5);
788 if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
789
790 Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
791 Int_t ip = gm->GetN();
792 gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
793 gm->SetPointError(ip, 0., 10.*f.GetParError(1));
794 gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
795 gs->SetPointError(ip, 0., 10.*f.GetParError(2));
796 }
797
798 // tracklet z resolution
799 h2 = (TH2I*)fContainer->At(kMCtrackletZ);
800 gm = (TGraphErrors*)fGraphM->At(kMCtrackletZ);
801 gs = (TGraphErrors*)fGraphS->At(kMCtrackletZ);
802 for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
803 h = h2->ProjectionY("py", iphi, iphi);
804 if(h->GetEntries()<100) continue;
805 AdjustF1(h, &fb);
806
807 if(IsVisual()){c->cd(); c->SetLogy();}
808 h->Fit(&fb, opt, "", -0.5, 0.5);
809 if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
810
811 Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
812 Int_t ip = gm->GetN();
813 gm->SetPoint(ip, phi, 10.*fb.GetParameter(1));
814 gm->SetPointError(ip, 0., 10.*fb.GetParError(1));
815 gs->SetPoint(ip, phi, 10.*fb.GetParameter(2));
816 gs->SetPointError(ip, 0., 10.*fb.GetParError(2));
817 }
818
819 //tracklet phi resolution
820 h2 = (TH2I*)fContainer->At(kMCtrackletPhi);
821 gm = (TGraphErrors*)fGraphM->At(kMCtrackletPhi);
822 gs = (TGraphErrors*)fGraphS->At(kMCtrackletPhi);
823 for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
824 h = h2->ProjectionY("py", iphi, iphi);
825 if(h->GetEntries()<100) continue;
826
827 if(IsVisual()){c->cd(); c->SetLogy();}
828 h->Fit(&f, opt, "", -0.5, 0.5);
829 if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
830
831 Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
832 Int_t ip = gm->GetN();
833 gm->SetPoint(ip, phi, f.GetParameter(1));
834 gm->SetPointError(ip, 0., f.GetParError(1));
835 gs->SetPoint(ip, phi, f.GetParameter(2));
836 gs->SetPointError(ip, 0., f.GetParError(2));
837 }
838
839 // track y resolution
840 h2 = (TH2I*)fContainer->At(kMCtrackY);
841 gm = (TGraphErrors*)fGraphM->At(kMCtrackY);
842 gs = (TGraphErrors*)fGraphS->At(kMCtrackY);
843 for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
844 h = h2->ProjectionY("py", iphi, iphi);
845 if(h->GetEntries()<100) continue;
846 AdjustF1(h, &f);
847
848 if(IsVisual()){c->cd(); c->SetLogy();}
849 h->Fit(&f, opt, "", -0.5, 0.5);
850 if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
851
852 Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
853 Int_t ip = gm->GetN();
854 gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
855 gm->SetPointError(ip, 0., 10.*f.GetParError(1));
856 gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
857 gs->SetPointError(ip, 0., 10.*f.GetParError(2));
d2381af5 858 }
765bd0ab 859
b1957d3c 860 // track z resolution
e3cf3d02 861 h2 = (TH2I*)fContainer->At(kMCtrackZIn);
862 gm = (TGraphErrors*)fGraphM->At(kMCtrackZIn);
863 gs = (TGraphErrors*)fGraphS->At(kMCtrackZIn);
b1957d3c 864 for(Int_t iphi=1; iphi<=h2->GetNbinsX(); iphi++){
865 h = h2->ProjectionY("pz", iphi, iphi);
866 if(h->GetEntries()<70) continue;
867 AdjustF1(h, &f);
868
869 if(IsVisual()){c->cd(); c->SetLogy();}
870 h->Fit(&f, opt, "", -0.5, 0.5);
871 if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
872
873 Double_t phi = h2->GetXaxis()->GetBinCenter(iphi);
874 Int_t ip = gm->GetN();
875 gm->SetPoint(ip, phi, 10.*f.GetParameter(1));
876 gm->SetPointError(ip, 0., 10.*f.GetParError(1));
877 gs->SetPoint(ip, phi, 10.*f.GetParameter(2));
878 gs->SetPointError(ip, 0., 10.*f.GetParError(2));
879 }
880
6090b78a 881 // track Pt resolution
882 h2 = (TH2I*)fContainer->At(kMCtrackPt);
883 TAxis *ax = h2->GetXaxis();
884 gm = (TGraphErrors*)fGraphM->At(kMCtrackPt);
885 gs = (TGraphErrors*)fGraphS->At(kMCtrackPt);
886 TF1 fg("fg", "gaus", -1.5, 1.5);
887 TF1 fl("fl", "landau", -4., 15.);
888 TF1 fgl("fgl", "gaus(0)+landau(3)", -5., 20.);
889 for(Int_t ip=1; ip<=ax->GetNbins(); ip++){
890 h = h2->ProjectionY("ppt", ip, ip);
891 if(h->GetEntries()<70) continue;
892
893 h->Fit(&fg, "QN", "", -1.5, 1.5);
894 fgl.SetParameter(0, fg.GetParameter(0));
895 fgl.SetParameter(1, fg.GetParameter(1));
896 fgl.SetParameter(2, fg.GetParameter(2));
897 h->Fit(&fl, "QN", "", -4., 15.);
898 fgl.SetParameter(3, fl.GetParameter(0));
899 fgl.SetParameter(4, fl.GetParameter(1));
900 fgl.SetParameter(5, fl.GetParameter(2));
901
902 if(IsVisual()){c->cd(); c->SetLogy();}
903 h->Fit(&fgl, opt, "", -5., 20.);
904 if(IsVisual()){c->Modified(); c->Update(); gSystem->Sleep(500);}
905
906 Float_t invpt = ax->GetBinCenter(ip);
907 Int_t ip = gm->GetN();
908 gm->SetPoint(ip, invpt, fgl.GetParameter(1));
909 gm->SetPointError(ip, 0., fgl.GetParError(1));
910 gs->SetPoint(ip, invpt, fgl.GetParameter(2)*invpt);
911 gs->SetPointError(ip, 0., fgl.GetParError(2));
912 // fgl.GetParameter(4) // Landau MPV
913 // fgl.GetParameter(5) // Landau Sigma
914 }
915
b1957d3c 916
917 if(c) delete c;
d85cd79c 918 return kTRUE;
919}
920
921
922//________________________________________________________
e3cf3d02 923void AliTRDresolution::Terminate(Option_t *)
d85cd79c 924{
925 if(fDebugStream){
926 delete fDebugStream;
927 fDebugStream = 0x0;
928 fDebugLevel = 0;
929 }
930 if(HasPostProcess()) PostProcess();
874acced 931}
d2381af5 932
3c3d9ff1 933//________________________________________________________
e3cf3d02 934void AliTRDresolution::AdjustF1(TH1 *h, TF1 *f)
3c3d9ff1 935{
936// Helper function to avoid duplication of code
937// Make first guesses on the fit parameters
938
939 // find the intial parameters of the fit !! (thanks George)
940 Int_t nbinsy = Int_t(.5*h->GetNbinsX());
941 Double_t sum = 0.;
942 for(Int_t jbin=nbinsy-4; jbin<=nbinsy+4; jbin++) sum+=h->GetBinContent(jbin); sum/=9.;
943 f->SetParLimits(0, 0., 3.*sum);
944 f->SetParameter(0, .9*sum);
945
017bd6af 946 f->SetParLimits(1, -.2, .2);
612ae7ed 947 f->SetParameter(1, -0.1);
3c3d9ff1 948
017bd6af 949 f->SetParLimits(2, 0., 4.e-1);
3c3d9ff1 950 f->SetParameter(2, 2.e-2);
017bd6af 951 if(f->GetNpar() <= 4) return;
3c3d9ff1 952
953 f->SetParLimits(3, 0., sum);
954 f->SetParameter(3, .1*sum);
955
956 f->SetParLimits(4, -.3, .3);
957 f->SetParameter(4, 0.);
958
959 f->SetParLimits(5, 0., 1.e2);
960 f->SetParameter(5, 2.e-1);
3c3d9ff1 961}
962
874acced 963//________________________________________________________
e3cf3d02 964TObjArray* AliTRDresolution::Histos()
874acced 965{
cf194b94 966 if(fContainer) return fContainer;
967
e3cf3d02 968 fContainer = new TObjArray(kNhistos);
94f34623 969 //fContainer->SetOwner(kTRUE);
cf194b94 970
b11eae29 971 TH1 *h = 0x0;
cf194b94 972 // cluster to tracklet residuals [2]
b1957d3c 973 if(!(h = (TH2I*)gROOT->FindObject("hCl"))){
974 h = new TH2I("hCl", "Clusters-Track Residuals", 21, -.33, .33, 100, -.5, .5);
124d488a 975 h->GetXaxis()->SetTitle("tg(#phi)");
976 h->GetYaxis()->SetTitle("#Delta y [cm]");
977 h->GetZaxis()->SetTitle("entries");
978 } else h->Reset();
b1957d3c 979 fContainer->AddAt(h, kCluster);
124d488a 980
251a1ae6 981 // tracklet to track residuals [2]
b1957d3c 982 if(!(h = (TH2I*)gROOT->FindObject("hTrkltY"))){
983 h = new TH2I("hTrkltY", "Tracklets-Track Residuals (Y)", 21, -.33, .33, 100, -.5, .5);
124d488a 984 h->GetXaxis()->SetTitle("tg(#phi)");
985 h->GetYaxis()->SetTitle("#Delta y [cm]");
986 h->GetZaxis()->SetTitle("entries");
987 } else h->Reset();
b1957d3c 988 fContainer->AddAt(h, kTrackletY);
124d488a 989
251a1ae6 990 // tracklet to track residuals angular [2]
b1957d3c 991 if(!(h = (TH2I*)gROOT->FindObject("hTrkltPhi"))){
992 h = new TH2I("hTrkltPhi", "Tracklets-Track Residuals (#Phi)", 21, -.33, .33, 100, -.5, .5);
124d488a 993 h->GetXaxis()->SetTitle("tg(#phi)");
994 h->GetYaxis()->SetTitle("#Delta phi [#circ]");
995 h->GetZaxis()->SetTitle("entries");
996 } else h->Reset();
b1957d3c 997 fContainer->AddAt(h, kTrackletPhi);
251a1ae6 998
cf194b94 999
1000 // Resolution histos
d937ad7a 1001 if(!HasMCdata()) return fContainer;
1002
1003 // cluster y resolution [0]
b1957d3c 1004 if(!(h = (TH2I*)gROOT->FindObject("hMCcl"))){
1005 h = new TH2I("hMCcl", "Cluster Resolution", 31, -.48, .48, 100, -.5, .5);
d937ad7a 1006 h->GetXaxis()->SetTitle("tg(#phi)");
1007 h->GetYaxis()->SetTitle("#Delta y [cm]");
1008 h->GetZaxis()->SetTitle("entries");
1009 } else h->Reset();
b1957d3c 1010 fContainer->AddAt(h, kMCcluster);
d937ad7a 1011
1012 // tracklet y resolution [0]
b1957d3c 1013 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltY"))){
1014 h = new TH2I("hMCtrkltY", "Tracklet Resolution (Y)", 31, -.48, .48, 100, -.5, .5);
d937ad7a 1015 h->GetXaxis()->SetTitle("tg(#phi)");
1016 h->GetYaxis()->SetTitle("#Delta y [cm]");
1017 h->GetZaxis()->SetTitle("entries");
1018 } else h->Reset();
b1957d3c 1019 fContainer->AddAt(h, kMCtrackletY);
d937ad7a 1020
e3cf3d02 1021 // tracklet y resolution [0]
1022 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltYPull"))){
1023 h = new TH2I("hMCtrkltYPull", "Tracklet Pulls (Y)", 31, -.48, .48, 100, -3.5, 3.5);
1024 h->GetXaxis()->SetTitle("tg(#phi)");
1025 h->GetYaxis()->SetTitle("#Delta y / #sigma_{y}");
1026 h->GetZaxis()->SetTitle("entries");
1027 } else h->Reset();
1028 fContainer->AddAt(h, kMCtrackletYPull);
1029
d937ad7a 1030 // tracklet y resolution [0]
b1957d3c 1031 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltZ"))){
1032 h = new TH2I("hMCtrkltZ", "Tracklet Resolution (Z)", 31, -.48, .48, 100, -.5, .5);
d937ad7a 1033 h->GetXaxis()->SetTitle("tg(#theta)");
1034 h->GetYaxis()->SetTitle("#Delta z [cm]");
1035 h->GetZaxis()->SetTitle("entries");
1036 } else h->Reset();
b1957d3c 1037 fContainer->AddAt(h, kMCtrackletZ);
d937ad7a 1038
e3cf3d02 1039 // tracklet y resolution [0]
1040 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltZPull"))){
1041 h = new TH2I("hMCtrkltZ", "Tracklet Pulls (Z)", 31, -.48, .48, 100, -3.5, 3.5);
1042 h->GetXaxis()->SetTitle("tg(#theta)");
1043 h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}");
1044 h->GetZaxis()->SetTitle("entries");
1045 } else h->Reset();
1046 fContainer->AddAt(h, kMCtrackletZPull);
1047
d937ad7a 1048 // tracklet angular resolution [1]
b1957d3c 1049 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkltPhi"))){
1050 h = new TH2I("hMCtrkltPhi", "Tracklet Resolution (#Phi)", 31, -.48, .48, 100, -10., 10.);
d937ad7a 1051 h->GetXaxis()->SetTitle("tg(#phi)");
1052 h->GetYaxis()->SetTitle("#Delta #phi [deg]");
1053 h->GetZaxis()->SetTitle("entries");
1054 } else h->Reset();
b1957d3c 1055 fContainer->AddAt(h, kMCtrackletPhi);
d937ad7a 1056
1057 // Kalman track y resolution
b1957d3c 1058 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkY"))){
1059 h = new TH2I("hMCtrkY", "Kalman Track Resolution (Y)", 31, -.48, .48, 100, -.5, .5);
d937ad7a 1060 h->GetXaxis()->SetTitle("tg(#phi)");
1061 h->GetYaxis()->SetTitle("#Delta y [cm]");
1062 h->GetZaxis()->SetTitle("entries");
1063 } else h->Reset();
b1957d3c 1064 fContainer->AddAt(h, kMCtrackY);
d937ad7a 1065
e3cf3d02 1066 // Kalman track y resolution
1067 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkYPull"))){
1068 h = new TH2I("hMCtrkYPull", "Kalman Track Pulls (Y)", 31, -.48, .48, 100, -3.5, 3.5);
1069 h->GetXaxis()->SetTitle("tg(#phi)");
1070 h->GetYaxis()->SetTitle("#Delta y / #sigma_{y}");
1071 h->GetZaxis()->SetTitle("entries");
1072 } else h->Reset();
1073 fContainer->AddAt(h, kMCtrackYPull);
1074
1075 // Kalman track Z resolution
1076 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZIn"))){
1077 h = new TH2I("hMCtrkZIn", "Kalman Track Resolution (Zin)", 20, -1., 1., 100, -1.5, 1.5);
1078 h->GetXaxis()->SetTitle("tg(#theta)");
1079 h->GetYaxis()->SetTitle("#Delta z [cm]");
1080 h->GetZaxis()->SetTitle("entries");
1081 } else h->Reset();
1082 fContainer->AddAt(h, kMCtrackZIn);
1083
d937ad7a 1084 // Kalman track Z resolution
e3cf3d02 1085 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZOut"))){
1086 h = new TH2I("hMCtrkZOut", "Kalman Track Resolution (Zout)", 20, -1., 1., 100, -1.5, 1.5);
d937ad7a 1087 h->GetXaxis()->SetTitle("tg(#theta)");
1088 h->GetYaxis()->SetTitle("#Delta z [cm]");
1089 h->GetZaxis()->SetTitle("entries");
1090 } else h->Reset();
e3cf3d02 1091 fContainer->AddAt(h, kMCtrackZOut);
1092
1093 // Kalman track Z resolution
1094 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZInPull"))){
1095 h = new TH2I("hMCtrkZInPull", "Kalman Track Pulls (Zin)", 20, -1., 1., 100, -4.5, 4.5);
1096 h->GetXaxis()->SetTitle("tg(#theta)");
1097 h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}");
1098 h->GetZaxis()->SetTitle("entries");
1099 } else h->Reset();
1100 fContainer->AddAt(h, kMCtrackZInPull);
1101
1102 // Kalman track Z resolution
1103 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkZOutPull"))){
1104 h = new TH2I("hMCtrkZOutPull", "Kalman Track Pulls (Zout)", 20, -1., 1., 100, -4.5, 4.5);
1105 h->GetXaxis()->SetTitle("tg(#theta)");
1106 h->GetYaxis()->SetTitle("#Delta z / #sigma_{z}");
1107 h->GetZaxis()->SetTitle("entries");
1108 } else h->Reset();
1109 fContainer->AddAt(h, kMCtrackZOutPull);
d937ad7a 1110
6090b78a 1111 // Kalman track Pt resolution
1112 if(!(h = (TH2I*)gROOT->FindObject("hMCtrkPt"))){
1113 h = new TH2I("hMCtrkPt", "Kalman Track Resolution (Pt)", 100, 0., 2., 150, -5., 20.);
1114 h->GetXaxis()->SetTitle("1/p_{t}");
1115 h->GetYaxis()->SetTitle("#Delta p_{t} [GeV/c]");
1116 h->GetZaxis()->SetTitle("entries");
1117 } else h->Reset();
1118 fContainer->AddAt(h, kMCtrackPt);
1119
3d86166d 1120 return fContainer;
77203477 1121}
1122
aaf47b30 1123
1124//________________________________________________________
e3cf3d02 1125void AliTRDresolution::SetRecoParam(AliTRDrecoParam *r)
aaf47b30 1126{
3d86166d 1127
aaf47b30 1128 fReconstructor->SetRecoParam(r);
1129}