]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG1/TRD/AliTRDresolution.cxx
Update to the analysis on MC data only.
[u/mrichter/AliRoot.git] / PWG1 / TRD / AliTRDresolution.cxx
CommitLineData
1ee39b3a 1/**************************************************************************
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**************************************************************************/
15
16/* $Id: AliTRDresolution.cxx 27496 2008-07-22 08:35:45Z cblume $ */
17
18////////////////////////////////////////////////////////////////////////////
19// //
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");
36// AliTRDresolution *res = new AliTRDresolution();
37// //res->SetMCdata();
38// //res->SetVerbose();
39// //res->SetVisual();
40// res->Load();
41// if(!res->PostProcess()) return;
42// res->GetRefFigure(0);
43// }
44//
45// Authors: //
46// Alexandru Bercuci <A.Bercuci@gsi.de> //
47// Markus Fasel <M.Fasel@gsi.de> //
48// //
49////////////////////////////////////////////////////////////////////////////
50
92c40e64 51#include <TSystem.h>
52
1ee39b3a 53#include <TROOT.h>
54#include <TObjArray.h>
55#include <TH3.h>
56#include <TH2.h>
57#include <TH1.h>
58#include <TF1.h>
59#include <TCanvas.h>
60#include <TGaxis.h>
61#include <TBox.h>
92c40e64 62#include <TLegend.h>
1ee39b3a 63#include <TGraphErrors.h>
64#include <TGraphAsymmErrors.h>
65#include <TMath.h>
66#include <TMatrixT.h>
67#include <TVectorT.h>
68#include <TTreeStream.h>
69#include <TGeoManager.h>
92c40e64 70#include <TDatabasePDG.h>
1ee39b3a 71
72#include "AliPID.h"
92d6d80c 73#include "AliLog.h"
92c40e64 74#include "AliESDtrack.h"
801d4d50 75#include "AliCDBManager.h"
76#include "AliCDBPath.h"
77#include "AliCDBEntry.h"
78#include "AliGeomManager.h"
068e2c00 79#include "AliMathBase.h"
1ee39b3a 80
81#include "AliTRDresolution.h"
82#include "AliTRDgeometry.h"
83#include "AliTRDpadPlane.h"
84#include "AliTRDcluster.h"
85#include "AliTRDseedV1.h"
86#include "AliTRDtrackV1.h"
87#include "AliTRDReconstructor.h"
88#include "AliTRDrecoParam.h"
92c40e64 89#include "AliTRDpidUtil.h"
1ee39b3a 90
91#include "info/AliTRDclusterInfo.h"
92
93ClassImp(AliTRDresolution)
94
92c40e64 95UChar_t const AliTRDresolution::fgNproj[kNviews] = {
a310e49b 96 2, 2, 5, 5, 5,
92d6d80c 97 2, 5, 11, 11, 11
1ee39b3a 98};
92c40e64 99Char_t const * AliTRDresolution::fgPerformanceName[kNviews] = {
1ee39b3a 100 "Charge"
101 ,"Cluster2Track"
102 ,"Tracklet2Track"
a310e49b 103 ,"Tracklet2TRDin"
104 ,"Tracklet2TRDout"
1ee39b3a 105 ,"Cluster2MC"
106 ,"Tracklet2MC"
92d6d80c 107 ,"TRDin2MC"
108 ,"TRDout2MC"
1ee39b3a 109 ,"TRD2MC"
110};
041d572e 111// Configure segmentation for y resolution/residuals
3ba3e0a4 112Int_t const AliTRDresolution::fgkNresYsegm[3] = {
113 AliTRDgeometry::kNsector
114 ,AliTRDgeometry::kNsector*AliTRDgeometry::kNstack
115 ,AliTRDgeometry::kNdet
116};
117Char_t const *AliTRDresolution::fgkResYsegmName[3] = {
118 "Sector", "Stack", "Detector"};
041d572e 119
1ee39b3a 120
121//________________________________________________________
122AliTRDresolution::AliTRDresolution()
f2e89a4c 123 :AliTRDrecoTask()
2589cf64 124 ,fSegmentLevel(0)
1ee39b3a 125 ,fIdxPlot(0)
92c40e64 126 ,fIdxFrame(0)
3ba3e0a4 127 ,fPtThreshold(1.)
b91fdd71 128 ,fReconstructor(NULL)
129 ,fGeo(NULL)
92c40e64 130 ,fDBPDG(NULL)
b91fdd71 131 ,fGraphS(NULL)
132 ,fGraphM(NULL)
133 ,fCl(NULL)
b91fdd71 134 ,fMCcl(NULL)
83b44483 135/* ,fTrklt(NULL)
136 ,fMCtrklt(NULL)*/
f8f46e4d 137{
138 //
139 // Default constructor
140 //
f2e89a4c 141 SetNameTitle("TRDresolution", "TRD spatial and momentum resolution");
2589cf64 142 SetSegmentationLevel();
f8f46e4d 143}
144
705f8b0a 145//________________________________________________________
f8f46e4d 146AliTRDresolution::AliTRDresolution(char* name)
f2e89a4c 147 :AliTRDrecoTask(name, "TRD spatial and momentum resolution")
2589cf64 148 ,fSegmentLevel(0)
f8f46e4d 149 ,fIdxPlot(0)
150 ,fIdxFrame(0)
3ba3e0a4 151 ,fPtThreshold(1.)
f8f46e4d 152 ,fReconstructor(NULL)
153 ,fGeo(NULL)
154 ,fDBPDG(NULL)
155 ,fGraphS(NULL)
156 ,fGraphM(NULL)
157 ,fCl(NULL)
f8f46e4d 158 ,fMCcl(NULL)
83b44483 159/* ,fTrklt(NULL)
160 ,fMCtrklt(NULL)*/
1ee39b3a 161{
162 //
163 // Default constructor
164 //
165
166 fReconstructor = new AliTRDReconstructor();
167 fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
168 fGeo = new AliTRDgeometry();
169
170 InitFunctorList();
2589cf64 171 SetSegmentationLevel();
1ee39b3a 172
705f8b0a 173 DefineOutput(kClToTrk, TObjArray::Class()); // cluster2track
705f8b0a 174 DefineOutput(kClToMC, TObjArray::Class()); // cluster2mc
83b44483 175/* DefineOutput(kTrkltToTrk, TObjArray::Class()); // tracklet2track
176 DefineOutput(kTrkltToMC, TObjArray::Class()); // tracklet2mc*/
1ee39b3a 177}
178
179//________________________________________________________
180AliTRDresolution::~AliTRDresolution()
181{
182 //
183 // Destructor
184 //
185
186 if(fGraphS){fGraphS->Delete(); delete fGraphS;}
187 if(fGraphM){fGraphM->Delete(); delete fGraphM;}
188 delete fGeo;
189 delete fReconstructor;
190 if(gGeoManager) delete gGeoManager;
191 if(fCl){fCl->Delete(); delete fCl;}
1ee39b3a 192 if(fMCcl){fMCcl->Delete(); delete fMCcl;}
83b44483 193/* if(fTrklt){fTrklt->Delete(); delete fTrklt;}
194 if(fMCtrklt){fMCtrklt->Delete(); delete fMCtrklt;}*/
1ee39b3a 195}
196
197
198//________________________________________________________
f8f46e4d 199void AliTRDresolution::UserCreateOutputObjects()
1ee39b3a 200{
201 // spatial resolution
5935a6da 202
fe1d1beb 203 if(!fReconstructor){
204 fReconstructor = new AliTRDReconstructor();
205 fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
206 }
fe1d1beb 207
068e2c00 208 AliTRDrecoTask::UserCreateOutputObjects();
209 PostData(kClToTrk, fCl);
210 PostData(kClToMC, fMCcl);
211/* PostData(kTrkltToTrk, fTrklt);
212 PostData(kTrkltToMC, fMCtrklt);*/
83b44483 213 InitExchangeContainers();
214}
1ee39b3a 215
83b44483 216//________________________________________________________
217void AliTRDresolution::InitExchangeContainers()
218{
1ee39b3a 219 fCl = new TObjArray();
220 fCl->SetOwner(kTRUE);
1ee39b3a 221 fMCcl = new TObjArray();
222 fMCcl->SetOwner(kTRUE);
83b44483 223/* fTrklt = new TObjArray();
224 fTrklt->SetOwner(kTRUE);
1ee39b3a 225 fMCtrklt = new TObjArray();
83b44483 226 fMCtrklt->SetOwner(kTRUE);*/
1ee39b3a 227}
228
229//________________________________________________________
f8f46e4d 230void AliTRDresolution::UserExec(Option_t *opt)
1ee39b3a 231{
232 //
233 // Execution part
234 //
235
801d4d50 236 if(!IsInitGeom()){
237 // create geometry
238 AliCDBManager* ocdb = AliCDBManager::Instance();
239 AliCDBEntry* obj = ocdb->Get(AliCDBPath("GRP", "Geometry", "Data"));
240 AliGeomManager::SetGeometry((TGeoManager*)obj->GetObject());
241 AliGeomManager::GetNalignable("TRD");
242 AliGeomManager::ApplyAlignObjsFromCDB("TRD");
243 fGeo = new AliTRDgeometry();
244 SetInitGeom();
245 }
246
1ee39b3a 247 fCl->Delete();
1ee39b3a 248 fMCcl->Delete();
83b44483 249/* fTrklt->Delete();
250 fMCtrklt->Delete();*/
b4414720 251 AliTRDrecoTask::UserExec(opt);
1ee39b3a 252}
253
553528eb 254//________________________________________________________
255Bool_t AliTRDresolution::Pulls(Double_t dyz[2], Double_t cov[3], Double_t tilt)
256{
257// Helper function to calculate pulls in the yz plane
258// using proper tilt rotation
259// Uses functionality defined by AliTRDseedV1.
260
261 Double_t t2(tilt*tilt);
262
263 // rotate along pad
264 Double_t cc[3];
265 cc[0] = cov[0] - 2.*tilt*cov[1] + t2*cov[2];
266 cc[1] = cov[1]*(1.-t2) + tilt*(cov[0] - cov[2]);
267 cc[2] = t2*cov[0] + 2.*tilt*cov[1] + cov[2];
268 // do sqrt
269 Double_t sqr[3]={0., 0., 0.};
270 if(AliTRDseedV1::GetCovSqrt(cc, sqr)) return kFALSE;
271 Double_t invsqr[3]={0., 0., 0.};
272 if(AliTRDseedV1::GetCovInv(sqr, invsqr)<1.e-40) return kFALSE;
273 Double_t tmp(dyz[0]);
274 dyz[0] = invsqr[0]*tmp + invsqr[1]*dyz[1];
275 dyz[1] = invsqr[1]*tmp + invsqr[2]*dyz[1];
276 return kTRUE;
277}
278
1ee39b3a 279//________________________________________________________
280TH1* AliTRDresolution::PlotCharge(const AliTRDtrackV1 *track)
281{
282 //
283 // Plots the charge distribution
284 //
285
286 if(track) fkTrack = track;
287 if(!fkTrack){
3d2a3dff 288 AliDebug(4, "No Track defined.");
b91fdd71 289 return NULL;
1ee39b3a 290 }
b91fdd71 291 TObjArray *arr = NULL;
83b44483 292 if(!fContainer || !(arr = ((TObjArray*)fContainer->At(kCharge)))){
1ee39b3a 293 AliWarning("No output container defined.");
b91fdd71 294 return NULL;
1ee39b3a 295 }
b91fdd71 296 TH3S* h = NULL;
1ee39b3a 297
b91fdd71 298 AliTRDseedV1 *fTracklet = NULL;
299 AliTRDcluster *c = NULL;
1ee39b3a 300 for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
301 if(!(fTracklet = fkTrack->GetTracklet(ily))) continue;
302 if(!fTracklet->IsOK()) continue;
303 Float_t x0 = fTracklet->GetX0();
304 Float_t dq, dl;
305 for(Int_t itb=AliTRDseedV1::kNtb; itb--;){
306 if(!(c = fTracklet->GetClusters(itb))){
307 if(!(c = fTracklet->GetClusters(AliTRDseedV1::kNtb+itb))) continue;
308 }
309 dq = fTracklet->GetdQdl(itb, &dl);
310 dl /= 0.15; // dl/dl0, dl0 = 1.5 mm for nominal vd
311 (h = (TH3S*)arr->At(0))->Fill(dl, x0-c->GetX(), dq);
312 }
313
314// if(!HasMCdata()) continue;
315// UChar_t s;
316// Float_t pt0, y0, z0, dydx0, dzdx0;
317// if(!fMC->GetDirections(x0, y0, z0, dydx0, dzdx0, pt0, s)) continue;
318
319 }
320 return h;
321}
322
323
324//________________________________________________________
325TH1* AliTRDresolution::PlotCluster(const AliTRDtrackV1 *track)
326{
327 //
328 // Plot the cluster distributions
329 //
330
331 if(track) fkTrack = track;
332 if(!fkTrack){
3d2a3dff 333 AliDebug(4, "No Track defined.");
b91fdd71 334 return NULL;
1ee39b3a 335 }
b91fdd71 336 TObjArray *arr = NULL;
83b44483 337 if(!fContainer || !(arr = ((TObjArray*)fContainer->At(kCluster)))){
1ee39b3a 338 AliWarning("No output container defined.");
b91fdd71 339 return NULL;
1ee39b3a 340 }
05039659 341 ULong_t status = fkESD ? fkESD->GetStatus():0;
1ee39b3a 342
3ba3e0a4 343 Int_t sgm[3];
553528eb 344 Double_t covR[7], cov[3], dy[2], dz[2];
345 Float_t pt, x0, y0, z0, dydx, dzdx;
dfd7d48b 346 AliTRDseedV1 *fTracklet(NULL); TObjArray *clInfoArr(NULL);
1ee39b3a 347 for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
348 if(!(fTracklet = fkTrack->GetTracklet(ily))) continue;
349 if(!fTracklet->IsOK()) continue;
350 x0 = fTracklet->GetX0();
3d2a3dff 351 pt = fTracklet->GetPt();
3ba3e0a4 352 sgm[2] = fTracklet->GetDetector();
353 sgm[0] = AliTRDgeometry::GetSector(sgm[2]);
354 sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]);
355
1ee39b3a 356 // retrive the track angle with the chamber
357 y0 = fTracklet->GetYref(0);
358 z0 = fTracklet->GetZref(0);
359 dydx = fTracklet->GetYref(1);
360 dzdx = fTracklet->GetZref(1);
31c8fa8a 361 fTracklet->GetCovRef(covR);
a47a87c5 362 Double_t tilt(fTracklet->GetTilt())
363 ,t2(tilt*tilt)
364 ,corr(1./(1. + t2))
6558fd69 365 ,cost(TMath::Sqrt(corr));
b91fdd71 366 AliTRDcluster *c = NULL;
1ee39b3a 367 fTracklet->ResetClusterIter(kFALSE);
368 while((c = fTracklet->PrevCluster())){
369 Float_t xc = c->GetX();
370 Float_t yc = c->GetY();
371 Float_t zc = c->GetZ();
372 Float_t dx = x0 - xc;
373 Float_t yt = y0 - dx*dydx;
374 Float_t zt = z0 - dx*dzdx;
a47a87c5 375 dy[0] = yc-yt; dz[0]= zc-zt;
376
553528eb 377 // rotate along pad
a47a87c5 378 dy[1] = cost*(dy[0] - dz[0]*tilt);
379 dz[1] = cost*(dz[0] + dy[0]*tilt);
83b44483 380 if(pt>fPtThreshold && c->IsInChamber()) ((TH3S*)arr->At(0))->Fill(dydx, dy[1], sgm[fSegmentLevel]);
31c8fa8a 381
553528eb 382 // tilt rotation of covariance for clusters
383 Double_t sy2(c->GetSigmaY2()), sz2(c->GetSigmaZ2());
a47a87c5 384 cov[0] = (sy2+t2*sz2)*corr;
385 cov[1] = tilt*(sz2 - sy2)*corr;
386 cov[2] = (t2*sy2+sz2)*corr;
553528eb 387 // sum with track covariance
a47a87c5 388 cov[0]+=covR[0]; cov[1]+=covR[1]; cov[2]+=covR[2];
553528eb 389 Double_t dyz[2]= {dy[1], dz[1]};
390 Pulls(dyz, cov, tilt);
2589cf64 391 ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]);
1ee39b3a 392
dfd7d48b 393 // Get z-position with respect to anode wire
394 Int_t istk = fGeo->GetStack(c->GetDetector());
395 AliTRDpadPlane *pp = fGeo->GetPadPlane(ily, istk);
396 Float_t row0 = pp->GetRow0();
397 Float_t d = row0 - zt + pp->GetAnodeWireOffset();
398 d -= ((Int_t)(2 * d)) / 2.0;
399 if (d > 0.25) d = 0.5 - d;
400
401 AliTRDclusterInfo *clInfo = new AliTRDclusterInfo;
dfd7d48b 402 clInfo->SetCluster(c);
403 Float_t covF[] = {cov[0], cov[1], cov[2]};
404 clInfo->SetGlobalPosition(yt, zt, dydx, dzdx, covF);
405 clInfo->SetResolution(dy[1]);
406 clInfo->SetAnisochronity(d);
407 clInfo->SetDriftLength(dx);
408 clInfo->SetTilt(tilt);
83b44483 409 if(fCl) fCl->Add(clInfo);
410 else AliDebug(1, "Cl exchange container missing. Activate by calling \"InitExchangeContainers()\"");
411
3ba3e0a4 412 if(DebugLevel()>=1){
d14a8ac2 413 if(!clInfoArr){
414 clInfoArr=new TObjArray(AliTRDseedV1::kNclusters);
415 clInfoArr->SetOwner(kFALSE);
416 }
dfd7d48b 417 clInfoArr->Add(clInfo);
df0514f6 418 }
1ee39b3a 419 }
3ba3e0a4 420 if(DebugLevel()>=1 && clInfoArr){
dfd7d48b 421 (*DebugStream()) << "cluster"
422 <<"status=" << status
423 <<"clInfo.=" << clInfoArr
424 << "\n";
d14a8ac2 425 clInfoArr->Clear();
dfd7d48b 426 }
1ee39b3a 427 }
d14a8ac2 428 if(clInfoArr) delete clInfoArr;
81979445 429 return (TH3S*)arr->At(0);
1ee39b3a 430}
431
432
433//________________________________________________________
434TH1* AliTRDresolution::PlotTracklet(const AliTRDtrackV1 *track)
435{
436// Plot normalized residuals for tracklets to track.
437//
438// We start from the result that if X=N(|m|, |Cov|)
439// BEGIN_LATEX
440// (Cov^{-1})^{1/2}X = N((Cov^{-1})^{1/2}*|m|, |1|)
441// END_LATEX
442// in our case X=(y_trklt - y_trk z_trklt - z_trk) and |Cov| = |Cov_trklt| + |Cov_trk| at the radial
443// reference position.
444 if(track) fkTrack = track;
445 if(!fkTrack){
3d2a3dff 446 AliDebug(4, "No Track defined.");
b91fdd71 447 return NULL;
1ee39b3a 448 }
b91fdd71 449 TObjArray *arr = NULL;
83b44483 450 if(!fContainer || !(arr = (TObjArray*)fContainer->At(kTrack ))){
1ee39b3a 451 AliWarning("No output container defined.");
b91fdd71 452 return NULL;
1ee39b3a 453 }
454
3ba3e0a4 455 Int_t sgm[3];
1ee39b3a 456 Double_t cov[3], covR[7]/*, sqr[3], inv[3]*/;
3ba3e0a4 457 Double_t pt, phi, tht, x, dx, dy[2], dz[2];
a47a87c5 458 AliTRDseedV1 *fTracklet(NULL);
3ba3e0a4 459 for(Int_t il(0); il<AliTRDgeometry::kNlayer; il++){
1ee39b3a 460 if(!(fTracklet = fkTrack->GetTracklet(il))) continue;
461 if(!fTracklet->IsOK()) continue;
3ba3e0a4 462 sgm[2] = fTracklet->GetDetector();
463 sgm[0] = AliTRDgeometry::GetSector(sgm[2]);
464 sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]);
a47a87c5 465 x = fTracklet->GetX();
466 dx = fTracklet->GetX0() - x;
467 pt = fTracklet->GetPt();
3ba3e0a4 468 phi = fTracklet->GetYref(1);
469 tht = fTracklet->GetZref(1);
a47a87c5 470 // compute dy and dz
471 dy[0]= fTracklet->GetYref(0)-dx*fTracklet->GetYref(1) - fTracklet->GetY();
472 dz[0]= fTracklet->GetZref(0)-dx*fTracklet->GetZref(1) - fTracklet->GetZ();
473 Double_t tilt(fTracklet->GetTilt())
474 ,t2(tilt*tilt)
475 ,corr(1./(1. + t2))
6558fd69 476 ,cost(TMath::Sqrt(corr));
477 Bool_t rc(fTracklet->IsRowCross());
a47a87c5 478
479 // calculate residuals using tilt rotation
480 dy[1]= cost*(dy[0] - dz[0]*tilt);
481 dz[1]= cost*(dz[0] + dy[0]*tilt);
2589cf64 482 ((TH3S*)arr->At(0))->Fill(phi, dy[1], sgm[fSegmentLevel]+rc*fgkNresYsegm[fSegmentLevel]);
3ba3e0a4 483 ((TH3S*)arr->At(2))->Fill(tht, dz[1], rc);
a47a87c5 484
1ee39b3a 485 // compute covariance matrix
486 fTracklet->GetCovAt(x, cov);
487 fTracklet->GetCovRef(covR);
488 cov[0] += covR[0]; cov[1] += covR[1]; cov[2] += covR[2];
553528eb 489 Double_t dyz[2]= {dy[1], dz[1]};
81979445 490 Pulls(dyz, cov, tilt);
2589cf64 491 ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]);
3ba3e0a4 492 ((TH3S*)arr->At(3))->Fill(tht, dyz[1], rc);
a47a87c5 493
3ba3e0a4 494 Double_t dphi((phi-fTracklet->GetYfit(1))/(1-phi*fTracklet->GetYfit(1)));
495 Double_t dtht((tht-fTracklet->GetZfit(1))/(1-tht*fTracklet->GetZfit(1)));
496 ((TH2I*)arr->At(4))->Fill(phi, TMath::ATan(dphi));
dfd7d48b 497
498 if(DebugLevel()>=1){
499 UChar_t err(fTracklet->GetErrorMsg());
500 (*DebugStream()) << "tracklet"
3ba3e0a4 501 <<"pt=" << pt
502 <<"phi=" << phi
503 <<"tht=" << tht
504 <<"det=" << sgm[2]
505 <<"dy0=" << dy[0]
506 <<"dz0=" << dz[0]
dfd7d48b 507 <<"dy=" << dy[1]
508 <<"dz=" << dz[1]
509 <<"dphi="<< dphi
510 <<"dtht="<< dtht
511 <<"dyp=" << dyz[0]
512 <<"dzp=" << dyz[1]
513 <<"rc=" << rc
514 <<"err=" << err
515 << "\n";
516 }
1ee39b3a 517 }
518
519
520 return (TH2I*)arr->At(0);
521}
522
523
524//________________________________________________________
a310e49b 525TH1* AliTRDresolution::PlotTrackIn(const AliTRDtrackV1 *track)
1ee39b3a 526{
527// Store resolution/pulls of Kalman before updating with the TRD information
528// at the radial position of the first tracklet. The following points are used
529// for comparison
530// - the (y,z,snp) of the first TRD tracklet
531// - the (y, z, snp, tgl, pt) of the MC track reference
532//
533// Additionally the momentum resolution/pulls are calculated for usage in the
534// PID calculation.
535
536 if(track) fkTrack = track;
537 if(!fkTrack){
3d2a3dff 538 AliDebug(4, "No Track defined.");
b91fdd71 539 return NULL;
1ee39b3a 540 }
83b44483 541 TObjArray *arr = NULL;
542 if(!fContainer || !(arr = (TObjArray*)fContainer->At(kTrackIn))){
543 AliWarning("No output container defined.");
544 return NULL;
545 }
b91fdd71 546 AliExternalTrackParam *tin = NULL;
a310e49b 547 if(!(tin = fkTrack->GetTrackIn())){
1ee39b3a 548 AliWarning("Track did not entered TRD fiducial volume.");
b91fdd71 549 return NULL;
1ee39b3a 550 }
b91fdd71 551 TH1 *h = NULL;
1ee39b3a 552
553 Double_t x = tin->GetX();
6558fd69 554 AliTRDseedV1 *fTracklet = NULL;
1ee39b3a 555 for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
6558fd69 556 if(!(fTracklet = fkTrack->GetTracklet(ily))) continue;
1ee39b3a 557 break;
558 }
6558fd69 559 if(!fTracklet || TMath::Abs(x-fTracklet->GetX())>1.e-3){
a310e49b 560 AliWarning("Tracklet did not match Track.");
b91fdd71 561 return NULL;
1ee39b3a 562 }
3ba3e0a4 563 Int_t sgm[3];
564 sgm[2] = fTracklet->GetDetector();
565 sgm[0] = AliTRDgeometry::GetSector(sgm[2]);
566 sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]);
6558fd69 567 Double_t tilt(fTracklet->GetTilt())
568 ,t2(tilt*tilt)
569 ,corr(1./(1. + t2))
570 ,cost(TMath::Sqrt(corr));
81979445 571 Bool_t rc(fTracklet->IsRowCross());
6558fd69 572
1ee39b3a 573 const Int_t kNPAR(5);
574 Double_t parR[kNPAR]; memcpy(parR, tin->GetParameter(), kNPAR*sizeof(Double_t));
575 Double_t covR[3*kNPAR]; memcpy(covR, tin->GetCovariance(), 3*kNPAR*sizeof(Double_t));
6558fd69 576 Double_t cov[3]; fTracklet->GetCovAt(x, cov);
1ee39b3a 577
578 // define sum covariances
579 TMatrixDSym COV(kNPAR); TVectorD PAR(kNPAR);
580 Double_t *pc = &covR[0], *pp = &parR[0];
581 for(Int_t ir=0; ir<kNPAR; ir++, pp++){
582 PAR(ir) = (*pp);
583 for(Int_t ic = 0; ic<=ir; ic++,pc++){
584 COV(ir,ic) = (*pc); COV(ic,ir) = (*pc);
585 }
586 }
587 PAR[4] = TMath::Abs(PAR[4]); // remove sign of pt !!
588 //COV.Print(); PAR.Print();
589
590 //TODO Double_t dydx = TMath::Sqrt(1.-parR[2]*parR[2])/parR[2];
553528eb 591 Double_t dy[2]={parR[0] - fTracklet->GetY(), 0.}
592 ,dz[2]={parR[1] - fTracklet->GetZ(), 0.}
6558fd69 593 ,dphi(TMath::ASin(PAR[2])-TMath::ATan(fTracklet->GetYfit(1)));
594 // calculate residuals using tilt rotation
595 dy[1] = cost*(dy[0] - dz[0]*tilt);
596 dz[1] = cost*(dz[0] + dy[0]*tilt);
597
2589cf64 598 if(1./PAR[4]>fPtThreshold) ((TH3S*)arr->At(0))->Fill(fTracklet->GetYref(1), dy[1], sgm[fSegmentLevel]+rc*fgkNresYsegm[fSegmentLevel]);
81979445 599 ((TH3S*)arr->At(2))->Fill(fTracklet->GetZref(1), dz[1], rc);
6558fd69 600 ((TH2I*)arr->At(4))->Fill(fTracklet->GetYref(1), dphi);
601
553528eb 602 Double_t dyz[2] = {dy[1], dz[1]};
603 Double_t cc[3] = {COV(0,0)+cov[0], COV(0,1)+cov[1], COV(1,1)+cov[2]};
604 Pulls(dyz, cc, tilt);
2589cf64 605 ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]);
81979445 606 ((TH3S*)arr->At(3))->Fill(fTracklet->GetZref(1), dyz[1], rc);
6558fd69 607
1ee39b3a 608
609
610 // register reference histo for mini-task
611 h = (TH2I*)arr->At(0);
612
3ba3e0a4 613 if(DebugLevel()>=2){
1ee39b3a 614 (*DebugStream()) << "trackIn"
615 << "x=" << x
616 << "P=" << &PAR
617 << "C=" << &COV
618 << "\n";
619
6558fd69 620 Double_t y = fTracklet->GetY();
621 Double_t z = fTracklet->GetZ();
1ee39b3a 622 (*DebugStream()) << "trackletIn"
623 << "y=" << y
624 << "z=" << z
625 << "Vy=" << cov[0]
626 << "Cyz=" << cov[1]
627 << "Vz=" << cov[2]
628 << "\n";
629 }
630
631
632 if(!HasMCdata()) return h;
633 UChar_t s;
6558fd69 634 Float_t dx, pt0, x0=fTracklet->GetX0(), y0, z0, dydx0, dzdx0;
1ee39b3a 635 if(!fkMC->GetDirections(x0, y0, z0, dydx0, dzdx0, pt0, s)) return h;
92c40e64 636 Int_t pdg = fkMC->GetPDG(),
637 sIdx(AliTRDpidUtil::Pdg2Pid(TMath::Abs(pdg))+1), // species index
638 sign(0);
639 if(!fDBPDG) fDBPDG=TDatabasePDG::Instance();
640 TParticlePDG *ppdg(fDBPDG->GetParticle(pdg));
641 if(ppdg) sign = ppdg->Charge() > 0. ? 1 : -1;
642
1ee39b3a 643 // translate to reference radial position
644 dx = x0 - x; y0 -= dx*dydx0; z0 -= dx*dzdx0;
645 Float_t norm = 1./TMath::Sqrt(1.+dydx0*dydx0); // 1/sqrt(1+tg^2(phi))
646 //Fill MC info
647 TVectorD PARMC(kNPAR);
648 PARMC[0]=y0; PARMC[1]=z0;
649 PARMC[2]=dydx0*norm; PARMC[3]=dzdx0*norm;
650 PARMC[4]=1./pt0;
651
652// TMatrixDSymEigen eigen(COV);
653// TVectorD evals = eigen.GetEigenValues();
654// TMatrixDSym evalsm(kNPAR);
655// for(Int_t ir=0; ir<kNPAR; ir++) for(Int_t ic=0; ic<kNPAR; ic++) evalsm(ir,ic) = (ir==ic ? evals(ir): 0.);
656// TMatrixD evecs = eigen.GetEigenVectors();
657// TMatrixD sqrcov(evecs, TMatrixD::kMult, TMatrixD(evalsm, TMatrixD::kMult, evecs.T()));
658
659 // fill histos
83b44483 660 if(!(arr = (TObjArray*)fContainer->At(kMCtrackIn))) {
661 AliWarning("No MC container defined.");
662 return h;
663 }
664
1ee39b3a 665 // y resolution/pulls
2589cf64 666 if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, PARMC[0]-PAR[0], sgm[fSegmentLevel]);
667 ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], (PARMC[0]-PAR[0])/TMath::Sqrt(COV(0,0)), (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1)));
1ee39b3a 668 // z resolution/pulls
81979445 669 ((TH3S*)arr->At(2))->Fill(dzdx0, PARMC[1]-PAR[1], 0);
670 ((TH3S*)arr->At(3))->Fill(dzdx0, (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1)), 0);
1ee39b3a 671 // phi resolution/snp pulls
672 ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ASin(PARMC[2])-TMath::ASin(PAR[2]));
673 ((TH2I*)arr->At(5))->Fill(dydx0, (PARMC[2]-PAR[2])/TMath::Sqrt(COV(2,2)));
674 // theta resolution/tgl pulls
675 ((TH2I*)arr->At(6))->Fill(dzdx0, TMath::ATan((PARMC[3]-PAR[3])/(1-PARMC[3]*PAR[3])));
676 ((TH2I*)arr->At(7))->Fill(dzdx0, (PARMC[3]-PAR[3])/TMath::Sqrt(COV(3,3)));
677 // pt resolution\\1/pt pulls\\p resolution/pull
92c40e64 678 ((TH3S*)arr->At(8))->Fill(pt0, PARMC[4]/PAR[4]-1., sign*sIdx);
679 ((TH3S*)arr->At(9))->Fill(PARMC[4], (PARMC[4]-PAR[4])/TMath::Sqrt(COV(4,4)), sign*sIdx);
680
681 Double_t p0 = TMath::Sqrt(1.+ PARMC[3]*PARMC[3])*pt0, p;
4226db3e 682 p = TMath::Sqrt(1.+ PAR[3]*PAR[3])/PAR[4];
92c40e64 683 ((TH3S*)arr->At(10))->Fill(p0, p/p0-1., sign*sIdx);
afca20ef 684// Float_t sp =
685// p*p*PAR[4]*PAR[4]*COV(4,4)
686// +2.*PAR[3]*COV(3,4)/PAR[4]
687// +PAR[3]*PAR[3]*COV(3,3)/p/p/PAR[4]/PAR[4]/PAR[4]/PAR[4];
688// if(sp>0.) ((TH3S*)arr->At(11))->Fill(p0, (p0-p)/TMath::Sqrt(sp), sign*sIdx);
1ee39b3a 689
690 // fill debug for MC
3ba3e0a4 691 if(DebugLevel()>=3){
1ee39b3a 692 (*DebugStream()) << "trackInMC"
693 << "P=" << &PARMC
694 << "\n";
695 }
696 return h;
697}
698
a310e49b 699//________________________________________________________
700TH1* AliTRDresolution::PlotTrackOut(const AliTRDtrackV1 *track)
701{
702// Store resolution/pulls of Kalman after last update with the TRD information
703// at the radial position of the first tracklet. The following points are used
704// for comparison
705// - the (y,z,snp) of the first TRD tracklet
706// - the (y, z, snp, tgl, pt) of the MC track reference
707//
708// Additionally the momentum resolution/pulls are calculated for usage in the
709// PID calculation.
710
711 if(track) fkTrack = track;
712 if(!fkTrack){
713 AliDebug(4, "No Track defined.");
714 return NULL;
715 }
83b44483 716 TObjArray *arr = NULL;
717 if(!fContainer || !(arr = (TObjArray*)fContainer->At(kTrackOut))){
718 AliWarning("No output container defined.");
719 return NULL;
720 }
a310e49b 721 AliExternalTrackParam *tout = NULL;
722 if(!(tout = fkTrack->GetTrackOut())){
e9d62bd2 723 AliDebug(2, "Track did not exit TRD.");
a310e49b 724 return NULL;
725 }
726 TH1 *h(NULL);
727
728 Double_t x = tout->GetX();
6558fd69 729 AliTRDseedV1 *fTracklet(NULL);
a310e49b 730 for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
6558fd69 731 if(!(fTracklet = fkTrack->GetTracklet(ily))) continue;
a310e49b 732 break;
733 }
6558fd69 734 if(!fTracklet || TMath::Abs(x-fTracklet->GetX())>1.e-3){
a310e49b 735 AliWarning("Tracklet did not match Track position.");
736 return NULL;
737 }
3ba3e0a4 738 Int_t sgm[3];
739 sgm[2] = fTracklet->GetDetector();
740 sgm[0] = AliTRDgeometry::GetSector(sgm[2]);
741 sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]);
6558fd69 742 Double_t tilt(fTracklet->GetTilt())
743 ,t2(tilt*tilt)
744 ,corr(1./(1. + t2))
745 ,cost(TMath::Sqrt(corr));
81979445 746 Bool_t rc(fTracklet->IsRowCross());
6558fd69 747
a310e49b 748 const Int_t kNPAR(5);
749 Double_t parR[kNPAR]; memcpy(parR, tout->GetParameter(), kNPAR*sizeof(Double_t));
750 Double_t covR[3*kNPAR]; memcpy(covR, tout->GetCovariance(), 3*kNPAR*sizeof(Double_t));
6558fd69 751 Double_t cov[3]; fTracklet->GetCovAt(x, cov);
a310e49b 752
753 // define sum covariances
754 TMatrixDSym COV(kNPAR); TVectorD PAR(kNPAR);
755 Double_t *pc = &covR[0], *pp = &parR[0];
756 for(Int_t ir=0; ir<kNPAR; ir++, pp++){
757 PAR(ir) = (*pp);
758 for(Int_t ic = 0; ic<=ir; ic++,pc++){
759 COV(ir,ic) = (*pc); COV(ic,ir) = (*pc);
760 }
761 }
762 PAR[4] = TMath::Abs(PAR[4]); // remove sign of pt !!
763 //COV.Print(); PAR.Print();
764
765 //TODO Double_t dydx = TMath::Sqrt(1.-parR[2]*parR[2])/parR[2];
6558fd69 766 Double_t dy[3]={parR[0] - fTracklet->GetY(), 0., 0.}
767 ,dz[3]={parR[1] - fTracklet->GetZ(), 0., 0.}
768 ,dphi(TMath::ASin(PAR[2])-TMath::ATan(fTracklet->GetYfit(1)));
769 // calculate residuals using tilt rotation
770 dy[1] = cost*(dy[0] - dz[0]*tilt);
771 dz[1] = cost*(dz[0] + dy[0]*tilt);
6558fd69 772
2589cf64 773 if(1./PAR[4]>fPtThreshold) ((TH3S*)arr->At(0))->Fill(fTracklet->GetYref(1), 1.e2*dy[1], sgm[fSegmentLevel]+rc*fgkNresYsegm[fSegmentLevel]); // scale to fit general residual range !!!
81979445 774 ((TH3S*)arr->At(2))->Fill(fTracklet->GetZref(1), dz[1], rc);
6558fd69 775 ((TH2I*)arr->At(4))->Fill(fTracklet->GetYref(1), dphi);
776
553528eb 777 Double_t dyz[2] = {dy[1], dz[1]};
778 Double_t cc[3] = {COV(0,0)+cov[0], COV(0,1)+cov[1], COV(1,1)+cov[2]};
779 Pulls(dyz, cc, tilt);
2589cf64 780 ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dyz[0], dyz[1]);
81979445 781 ((TH3S*)arr->At(3))->Fill(fTracklet->GetZref(1), dyz[1], rc);
a310e49b 782
783 // register reference histo for mini-task
784 h = (TH2I*)arr->At(0);
785
3ba3e0a4 786 if(DebugLevel()>=2){
a310e49b 787 (*DebugStream()) << "trackOut"
788 << "x=" << x
789 << "P=" << &PAR
790 << "C=" << &COV
791 << "\n";
792
6558fd69 793 Double_t y = fTracklet->GetY();
794 Double_t z = fTracklet->GetZ();
a310e49b 795 (*DebugStream()) << "trackletOut"
796 << "y=" << y
797 << "z=" << z
798 << "Vy=" << cov[0]
799 << "Cyz=" << cov[1]
800 << "Vz=" << cov[2]
801 << "\n";
802 }
803
804
805 if(!HasMCdata()) return h;
806 UChar_t s;
6558fd69 807 Float_t dx, pt0, x0=fTracklet->GetX0(), y0, z0, dydx0, dzdx0;
a310e49b 808 if(!fkMC->GetDirections(x0, y0, z0, dydx0, dzdx0, pt0, s)) return h;
809 Int_t pdg = fkMC->GetPDG(),
810 sIdx(AliTRDpidUtil::Pdg2Pid(TMath::Abs(pdg))+1), // species index
811 sign(0);
812 if(!fDBPDG) fDBPDG=TDatabasePDG::Instance();
813 TParticlePDG *ppdg(fDBPDG->GetParticle(pdg));
814 if(ppdg) sign = ppdg->Charge() > 0. ? 1 : -1;
815
816 // translate to reference radial position
817 dx = x0 - x; y0 -= dx*dydx0; z0 -= dx*dzdx0;
818 Float_t norm = 1./TMath::Sqrt(1.+dydx0*dydx0); // 1/sqrt(1+tg^2(phi))
819 //Fill MC info
820 TVectorD PARMC(kNPAR);
821 PARMC[0]=y0; PARMC[1]=z0;
822 PARMC[2]=dydx0*norm; PARMC[3]=dzdx0*norm;
823 PARMC[4]=1./pt0;
824
825// TMatrixDSymEigen eigen(COV);
826// TVectorD evals = eigen.GetEigenValues();
827// TMatrixDSym evalsm(kNPAR);
828// for(Int_t ir=0; ir<kNPAR; ir++) for(Int_t ic=0; ic<kNPAR; ic++) evalsm(ir,ic) = (ir==ic ? evals(ir): 0.);
829// TMatrixD evecs = eigen.GetEigenVectors();
830// TMatrixD sqrcov(evecs, TMatrixD::kMult, TMatrixD(evalsm, TMatrixD::kMult, evecs.T()));
831
832 // fill histos
83b44483 833 if(!(arr = (TObjArray*)fContainer->At(kMCtrackOut))){
834 AliWarning("No MC container defined.");
835 return h;
836 }
a310e49b 837 // y resolution/pulls
2589cf64 838 if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, PARMC[0]-PAR[0], sgm[fSegmentLevel]);
839 ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], (PARMC[0]-PAR[0])/TMath::Sqrt(COV(0,0)), (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1)));
a310e49b 840 // z resolution/pulls
81979445 841 ((TH3S*)arr->At(2))->Fill(dzdx0, PARMC[1]-PAR[1], 0);
842 ((TH3S*)arr->At(3))->Fill(dzdx0, (PARMC[1]-PAR[1])/TMath::Sqrt(COV(1,1)), 0);
a310e49b 843 // phi resolution/snp pulls
844 ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ASin(PARMC[2])-TMath::ASin(PAR[2]));
845 ((TH2I*)arr->At(5))->Fill(dydx0, (PARMC[2]-PAR[2])/TMath::Sqrt(COV(2,2)));
846 // theta resolution/tgl pulls
847 ((TH2I*)arr->At(6))->Fill(dzdx0, TMath::ATan((PARMC[3]-PAR[3])/(1-PARMC[3]*PAR[3])));
848 ((TH2I*)arr->At(7))->Fill(dzdx0, (PARMC[3]-PAR[3])/TMath::Sqrt(COV(3,3)));
849 // pt resolution\\1/pt pulls\\p resolution/pull
850 ((TH3S*)arr->At(8))->Fill(pt0, PARMC[4]/PAR[4]-1., sign*sIdx);
851 ((TH3S*)arr->At(9))->Fill(PARMC[4], (PARMC[4]-PAR[4])/TMath::Sqrt(COV(4,4)), sign*sIdx);
852
853 Double_t p0 = TMath::Sqrt(1.+ PARMC[3]*PARMC[3])*pt0, p;
854 p = TMath::Sqrt(1.+ PAR[3]*PAR[3])/PAR[4];
855 ((TH3S*)arr->At(10))->Fill(p0, p/p0-1., sign*sIdx);
856// Float_t sp =
857// p*p*PAR[4]*PAR[4]*COV(4,4)
858// +2.*PAR[3]*COV(3,4)/PAR[4]
859// +PAR[3]*PAR[3]*COV(3,3)/p/p/PAR[4]/PAR[4]/PAR[4]/PAR[4];
860// if(sp>0.) ((TH3S*)arr->At(11))->Fill(p0, (p0-p)/TMath::Sqrt(sp), sign*sIdx);
861
862 // fill debug for MC
3ba3e0a4 863 if(DebugLevel()>=3){
a310e49b 864 (*DebugStream()) << "trackOutMC"
865 << "P=" << &PARMC
866 << "\n";
867 }
868 return h;
869}
870
1ee39b3a 871//________________________________________________________
872TH1* AliTRDresolution::PlotMC(const AliTRDtrackV1 *track)
873{
874 //
875 // Plot MC distributions
876 //
877
878 if(!HasMCdata()){
e9d62bd2 879 AliDebug(2, "No MC defined. Results will not be available.");
b91fdd71 880 return NULL;
1ee39b3a 881 }
882 if(track) fkTrack = track;
883 if(!fkTrack){
3d2a3dff 884 AliDebug(4, "No Track defined.");
b91fdd71 885 return NULL;
1ee39b3a 886 }
83b44483 887 if(!fContainer){
888 AliWarning("No output container defined.");
889 return NULL;
890 }
92c40e64 891 // retriev track characteristics
892 Int_t pdg = fkMC->GetPDG(),
893 sIdx(AliTRDpidUtil::Pdg2Pid(TMath::Abs(pdg))+1), // species index
894 sign(0),
3ba3e0a4 895 sgm[3],
92c40e64 896 label(fkMC->GetLabel());
897 if(!fDBPDG) fDBPDG=TDatabasePDG::Instance();
898 TParticlePDG *ppdg(fDBPDG->GetParticle(pdg));
899 if(ppdg) sign = ppdg->Charge() > 0. ? 1 : -1;
92c40e64 900
901 TObjArray *arr(NULL);TH1 *h(NULL);
1ee39b3a 902 UChar_t s;
1ee39b3a 903 Double_t xAnode, x, y, z, pt, dydx, dzdx, dzdl;
904 Float_t pt0, x0, y0, z0, dx, dy, dz, dydx0, dzdx0;
905 Double_t covR[7]/*, cov[3]*/;
906
3ba3e0a4 907 if(DebugLevel()>=3){
5066aa9a 908 TVectorD dX(12), dY(12), dZ(12), Pt(12), dPt(12), cCOV(12*15);
909 fkMC->PropagateKalman(&dX, &dY, &dZ, &Pt, &dPt, &cCOV);
1ee39b3a 910 (*DebugStream()) << "MCkalman"
4226db3e 911 << "pdg=" << pdg
912 << "dx=" << &dX
913 << "dy=" << &dY
914 << "dz=" << &dZ
5066aa9a 915 << "pt=" << &Pt
4226db3e 916 << "dpt=" << &dPt
917 << "cov=" << &cCOV
1ee39b3a 918 << "\n";
919 }
920
921 AliTRDReconstructor rec;
d43e2ad1 922 AliTRDseedV1 *fTracklet(NULL); TObjArray *clInfoArr(NULL);
1ee39b3a 923 for(Int_t ily=0; ily<AliTRDgeometry::kNlayer; ily++){
924 if(!(fTracklet = fkTrack->GetTracklet(ily)))/* ||
925 !fTracklet->IsOK())*/ continue;
926
3ba3e0a4 927 sgm[2] = fTracklet->GetDetector();
928 sgm[0] = AliTRDgeometry::GetSector(sgm[2]);
929 sgm[1] = sgm[0] * AliTRDgeometry::kNstack + AliTRDgeometry::GetStack(sgm[2]);
6558fd69 930 Double_t tilt(fTracklet->GetTilt())
931 ,t2(tilt*tilt)
932 ,corr(1./(1. + t2))
933 ,cost(TMath::Sqrt(corr));
1ee39b3a 934 x0 = fTracklet->GetX0();
935 //radial shift with respect to the MC reference (radial position of the pad plane)
936 x= fTracklet->GetX();
81979445 937 Bool_t rc(fTracklet->IsRowCross());
1ee39b3a 938 if(!fkMC->GetDirections(x0, y0, z0, dydx0, dzdx0, pt0, s)) continue;
939 xAnode = fTracklet->GetX0();
940
941 // MC track position at reference radial position
942 dx = x0 - x;
3ba3e0a4 943 if(DebugLevel()>=4){
1ee39b3a 944 (*DebugStream()) << "MC"
3ba3e0a4 945 << "det=" << sgm[2]
1ee39b3a 946 << "pdg=" << pdg
92c40e64 947 << "sgn=" << sign
1ee39b3a 948 << "pt=" << pt0
949 << "x=" << x0
950 << "y=" << y0
951 << "z=" << z0
952 << "dydx=" << dydx0
953 << "dzdx=" << dzdx0
954 << "\n";
955 }
956 Float_t ymc = y0 - dx*dydx0;
957 Float_t zmc = z0 - dx*dzdx0;
958 //p = pt0*TMath::Sqrt(1.+dzdx0*dzdx0); // pt -> p
959
960 // Kalman position at reference radial position
961 dx = xAnode - x;
962 dydx = fTracklet->GetYref(1);
963 dzdx = fTracklet->GetZref(1);
964 dzdl = fTracklet->GetTgl();
965 y = fTracklet->GetYref(0) - dx*dydx;
966 dy = y - ymc;
967 z = fTracklet->GetZref(0) - dx*dzdx;
968 dz = z - zmc;
969 pt = TMath::Abs(fTracklet->GetPt());
970 fTracklet->GetCovRef(covR);
971
92d6d80c 972 arr = (TObjArray*)((TObjArray*)fContainer->At(kMCtrack))->At(ily);
1ee39b3a 973 // y resolution/pulls
2589cf64 974 if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, dy, sgm[fSegmentLevel]);
975 ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dy/TMath::Sqrt(covR[0]), dz/TMath::Sqrt(covR[2]));
1ee39b3a 976 // z resolution/pulls
81979445 977 ((TH3S*)arr->At(2))->Fill(dzdx0, dz, 0);
978 ((TH3S*)arr->At(3))->Fill(dzdx0, dz/TMath::Sqrt(covR[2]), 0);
1ee39b3a 979 // phi resolution/ snp pulls
980 Double_t dtgp = (dydx - dydx0)/(1.- dydx*dydx0);
981 ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ATan(dtgp));
982 Double_t dsnp = dydx/TMath::Sqrt(1.+dydx*dydx) - dydx0/TMath::Sqrt(1.+dydx0*dydx0);
983 ((TH2I*)arr->At(5))->Fill(dydx0, dsnp/TMath::Sqrt(covR[3]));
984 // theta resolution/ tgl pulls
985 Double_t dzdl0 = dzdx0/TMath::Sqrt(1.+dydx0*dydx0),
986 dtgl = (dzdl - dzdl0)/(1.- dzdl*dzdl0);
987 ((TH2I*)arr->At(6))->Fill(dzdl0,
988 TMath::ATan(dtgl));
989 ((TH2I*)arr->At(7))->Fill(dzdl0, (dzdl - dzdl0)/TMath::Sqrt(covR[4]));
990 // pt resolution \\ 1/pt pulls \\ p resolution for PID
92c40e64 991 Double_t p0 = TMath::Sqrt(1.+ dzdl0*dzdl0)*pt0,
992 p = TMath::Sqrt(1.+ dzdl*dzdl)*pt;
afca20ef 993 ((TH3S*)((TObjArray*)arr->At(8)))->Fill(pt0, pt/pt0-1., sign*sIdx);
994 ((TH3S*)((TObjArray*)arr->At(9)))->Fill(1./pt0, (1./pt-1./pt0)/TMath::Sqrt(covR[6]), sign*sIdx);
995 ((TH3S*)((TObjArray*)arr->At(10)))->Fill(p0, p/p0-1., sign*sIdx);
1ee39b3a 996
997 // Fill Debug stream for Kalman track
3ba3e0a4 998 if(DebugLevel()>=4){
1ee39b3a 999 (*DebugStream()) << "MCtrack"
1000 << "pt=" << pt
1001 << "x=" << x
1002 << "y=" << y
1003 << "z=" << z
1004 << "dydx=" << dydx
1005 << "dzdx=" << dzdx
1006 << "s2y=" << covR[0]
1007 << "s2z=" << covR[2]
1008 << "\n";
1009 }
1010
1011 // recalculate tracklet based on the MC info
1012 AliTRDseedV1 tt(*fTracklet);
1013 tt.SetZref(0, z0 - (x0-xAnode)*dzdx0);
1014 tt.SetZref(1, dzdx0);
1015 tt.SetReconstructor(&rec);
26e4e234 1016 tt.Fit(1);
1ee39b3a 1017 x= tt.GetX();y= tt.GetY();z= tt.GetZ();
1018 dydx = tt.GetYfit(1);
1019 dx = x0 - x;
1020 ymc = y0 - dx*dydx0;
1021 zmc = z0 - dx*dzdx0;
6558fd69 1022 dy = y-ymc;
1023 dz = z-zmc;
81979445 1024 Float_t dphi = (dydx - dydx0);
1025 dphi /= (1.- dydx*dydx0);
6558fd69 1026
1ee39b3a 1027 // add tracklet residuals for y and dydx
1028 arr = (TObjArray*)fContainer->At(kMCtracklet);
81979445 1029
2589cf64 1030 if(pt0>fPtThreshold) ((TH3S*)arr->At(0))->Fill(dydx0, dy, sgm[fSegmentLevel]);
1031 if(tt.GetS2Y()>0. && tt.GetS2Z()>0.) ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dy/TMath::Sqrt(tt.GetS2Y()), dz/TMath::Sqrt(tt.GetS2Z()));
81979445 1032 ((TH3S*)arr->At(2))->Fill(dzdl0, dz, rc);
1033 if(tt.GetS2Z()>0.) ((TH3S*)arr->At(3))->Fill(dzdl0, dz/TMath::Sqrt(tt.GetS2Z()), rc);
1034 ((TH2I*)arr->At(4))->Fill(dydx0, TMath::ATan(dphi));
1ee39b3a 1035
1036 // Fill Debug stream for tracklet
3ba3e0a4 1037 if(DebugLevel()>=4){
1ee39b3a 1038 Float_t s2y = tt.GetS2Y();
1039 Float_t s2z = tt.GetS2Z();
1040 (*DebugStream()) << "MCtracklet"
1041 << "rc=" << rc
1042 << "x=" << x
1043 << "y=" << y
1044 << "z=" << z
1045 << "dydx=" << dydx
1046 << "s2y=" << s2y
1047 << "s2z=" << s2z
1048 << "\n";
1049 }
1050
3ba3e0a4 1051 AliTRDpadPlane *pp = fGeo->GetPadPlane(ily, AliTRDgeometry::GetStack(sgm[2]));
1ee39b3a 1052 Float_t zr0 = pp->GetRow0() + pp->GetAnodeWireOffset();
1ee39b3a 1053 //Double_t exb = AliTRDCommonParam::Instance()->GetOmegaTau(1.5);
1054
1055 arr = (TObjArray*)fContainer->At(kMCcluster);
b91fdd71 1056 AliTRDcluster *c = NULL;
83b44483 1057 tt.ResetClusterIter(kFALSE);
1058 while((c = tt.PrevCluster())){
1ee39b3a 1059 Float_t q = TMath::Abs(c->GetQ());
1060 x = c->GetX(); y = c->GetY();z = c->GetZ();
1061 dx = x0 - x;
1062 ymc= y0 - dx*dydx0;
1063 zmc= z0 - dx*dzdx0;
6558fd69 1064 dy = cost*(y - ymc - tilt*(z-zmc));
1065 dz = cost*(z - zmc + tilt*(y-ymc));
1066
1ee39b3a 1067 // Fill Histograms
83b44483 1068 if(q>20. && q<250. && pt0>fPtThreshold && c->IsInChamber()){
2589cf64 1069 ((TH3S*)arr->At(0))->Fill(dydx0, dy, sgm[fSegmentLevel]);
1070 ((TH3S*)arr->At(1))->Fill(sgm[fSegmentLevel], dy/TMath::Sqrt(c->GetSigmaY2()), dz/TMath::Sqrt(c->GetSigmaZ2()));
1ee39b3a 1071 }
1072
1073 // Fill calibration container
1074 Float_t d = zr0 - zmc;
1075 d -= ((Int_t)(2 * d)) / 2.0;
1076 if (d > 0.25) d = 0.5 - d;
1077 AliTRDclusterInfo *clInfo = new AliTRDclusterInfo;
1ee39b3a 1078 clInfo->SetCluster(c);
1079 clInfo->SetMC(pdg, label);
1080 clInfo->SetGlobalPosition(ymc, zmc, dydx0, dzdx0);
1081 clInfo->SetResolution(dy);
1082 clInfo->SetAnisochronity(d);
2589cf64 1083 clInfo->SetDriftLength(dx);
1ee39b3a 1084 clInfo->SetTilt(tilt);
83b44483 1085 if(fMCcl) fMCcl->Add(clInfo);
1086 else AliDebug(1, "MCcl exchange container missing. Activate by calling \"InitExchangeContainers()\"");
3ba3e0a4 1087 if(DebugLevel()>=5){
d14a8ac2 1088 if(!clInfoArr){
1089 clInfoArr=new TObjArray(AliTRDseedV1::kNclusters);
1090 clInfoArr->SetOwner(kFALSE);
1091 }
d43e2ad1 1092 clInfoArr->Add(clInfo);
1ee39b3a 1093 }
1094 }
d43e2ad1 1095 // Fill Debug Tree
3ba3e0a4 1096 if(DebugLevel()>=5 && clInfoArr){
d43e2ad1 1097 (*DebugStream()) << "MCcluster"
1098 <<"clInfo.=" << clInfoArr
1099 << "\n";
92c40e64 1100 clInfoArr->Clear();
d43e2ad1 1101 }
1ee39b3a 1102 }
92c40e64 1103 if(clInfoArr) delete clInfoArr;
1ee39b3a 1104 return h;
1105}
1106
1107
1108//________________________________________________________
1109Bool_t AliTRDresolution::GetRefFigure(Int_t ifig)
1110{
1111 //
1112 // Get the reference figures
1113 //
1114
1115 Float_t xy[4] = {0., 0., 0., 0.};
1116 if(!gPad){
1117 AliWarning("Please provide a canvas to draw results.");
1118 return kFALSE;
1119 }
041d572e 1120 Int_t selection[100], n(0), selStart(0); //
a310e49b 1121 Int_t ly0(0), dly(5);
1122 //Int_t ly0(1), dly(2); // used for SA
92d6d80c 1123 TList *l(NULL); TVirtualPad *pad(NULL);
1124 TGraphErrors *g(NULL);TGraphAsymmErrors *ga(NULL);
1ee39b3a 1125 switch(ifig){
a310e49b 1126 case 0: // charge resolution
1ee39b3a 1127 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1128 ((TVirtualPad*)l->At(0))->cd();
92d6d80c 1129 ga=((TGraphAsymmErrors*)((TObjArray*)fGraphM->At(kCharge))->At(0));
1130 if(ga->GetN()) ga->Draw("apl");
1ee39b3a 1131 ((TVirtualPad*)l->At(1))->cd();
92d6d80c 1132 g = ((TGraphErrors*)((TObjArray*)fGraphS->At(kCharge))->At(0));
1133 if(g->GetN()) g->Draw("apl");
1ee39b3a 1134 break;
a310e49b 1135 case 1: // cluster2track residuals
1ee39b3a 1136 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
6558fd69 1137 xy[0] = -.3; xy[1] = -100.; xy[2] = .3; xy[3] = 1000.;
a310e49b 1138 pad = (TVirtualPad*)l->At(0); pad->cd();
1139 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
2589cf64 1140 selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1141 if(!GetGraphArray(xy, kCluster, 0, 1, n, selection)) break;
1142 pad=(TVirtualPad*)l->At(1); pad->cd();
1143 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
2589cf64 1144 selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1145 if(!GetGraphArray(xy, kCluster, 0, 1, n, selection)) break;
1146 return kTRUE;
1147 case 2: // cluster2track residuals
1148 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1149 xy[0] = -.3; xy[1] = -100.; xy[2] = .3; xy[3] = 1000.;
1150 pad = (TVirtualPad*)l->At(0); pad->cd();
1151 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
2589cf64 1152 selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1153 if(!GetGraphArray(xy, kCluster, 0, 1, n, selection)) break;
2589cf64 1154 xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-0.5; xy[3] = 2.5;
a310e49b 1155 pad=(TVirtualPad*)l->At(1); pad->cd();
1156 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
81979445 1157 if(!GetGraphArray(xy, kCluster, 1, 1)) break;
1ee39b3a 1158 return kTRUE;
6558fd69 1159 case 3: // kTrack y
dfd7d48b 1160 gPad->Divide(3, 2, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
3ba3e0a4 1161 xy[0] = -.3; xy[1] = -20.; xy[2] = .3; xy[3] = 100.;
1ee39b3a 1162 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1163 selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1164 if(!GetGraphArray(xy, kTrack, 0, 1, n, selection)) break;
dfd7d48b 1165
6558fd69 1166 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1167 selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1168 if(!GetGraphArray(xy, kTrack, 0, 1, n, selection)) break;
dfd7d48b 1169
1170 ((TVirtualPad*)l->At(2))->cd();
2589cf64 1171 selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
dfd7d48b 1172 if(!GetGraphArray(xy, kTrack, 0, 1, n, selection)) break;
1173
1174 ((TVirtualPad*)l->At(3))->cd();
2589cf64 1175 selStart=fgkNresYsegm[fSegmentLevel]; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
dfd7d48b 1176 if(!GetGraphArray(xy, kTrack, 0, 1, n, selection, "[RC]")) break;
1177
1178 ((TVirtualPad*)l->At(4))->cd();
2589cf64 1179 selStart=fgkNresYsegm[fSegmentLevel]/3+fgkNresYsegm[fSegmentLevel]; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
dfd7d48b 1180 if(!GetGraphArray(xy, kTrack, 0, 1, n, selection, "[RC]")) break;
1181
1182 ((TVirtualPad*)l->At(5))->cd();
2589cf64 1183 selStart=2*fgkNresYsegm[fSegmentLevel]/3+fgkNresYsegm[fSegmentLevel]; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
dfd7d48b 1184 if(!GetGraphArray(xy, kTrack, 0, 1, n, selection, "[RC]")) break;
6558fd69 1185 return kTRUE;
dfd7d48b 1186 case 4: // kTrack z
6558fd69 1187 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
dfd7d48b 1188
1189 xy[0] = -1.; xy[1] = -150.; xy[2] = 1.; xy[3] = 1000.;
6558fd69 1190 ((TVirtualPad*)l->At(0))->cd();
dfd7d48b 1191 selection[0]=1;
1192 if(!GetGraphArray(xy, kTrack, 2, 1, 1, selection)) break;
1193
1194 xy[0] = -1.; xy[1] = -1500.; xy[2] = 1.; xy[3] = 10000.;
1ee39b3a 1195 ((TVirtualPad*)l->At(1))->cd();
dfd7d48b 1196 selection[0]=0;
1197 if(!GetGraphArray(xy, kTrack, 2, 1, 1, selection)) break;
1198
1ee39b3a 1199 return kTRUE;
dfd7d48b 1200 case 5: // kTrack pulls
1ee39b3a 1201 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
dfd7d48b 1202
2589cf64 1203 xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 2.5;
1ee39b3a 1204 ((TVirtualPad*)l->At(0))->cd();
dfd7d48b 1205 if(!GetGraphArray(xy, kTrack, 1, 1)) break;
1206
1ee39b3a 1207 xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
1208 ((TVirtualPad*)l->At(1))->cd();
81979445 1209 if(!GetGraphArray(xy, kTrack, 3, 1)) break;
1ee39b3a 1210 return kTRUE;
6558fd69 1211 case 6: // kTrack phi
1ee39b3a 1212 xy[0] = -.3; xy[1] = -5.; xy[2] = .3; xy[3] = 50.;
a310e49b 1213 if(GetGraph(&xy[0], kTrack , 4)) return kTRUE;
1ee39b3a 1214 break;
6558fd69 1215 case 7: // kTrackIn y
1ee39b3a 1216 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
a310e49b 1217 xy[0] = -.3; xy[1] = -1500.; xy[2] = .3; xy[3] = 5000.;
1ee39b3a 1218 pad = ((TVirtualPad*)l->At(0)); pad->cd();
1219 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
2589cf64 1220 selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1221 if(!GetGraphArray(xy, kTrackIn, 0, 1, n, selection)) break;
1222 pad=((TVirtualPad*)l->At(1)); pad->cd();
1223 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
2589cf64 1224 selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1225 if(!GetGraphArray(xy, kTrackIn, 0, 1, n, selection)) break;
1226 return kTRUE;
1227 case 8: // kTrackIn y
1228 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1229 xy[0] = -.3; xy[1] = -1500.; xy[2] = .3; xy[3] = 5000.;
1230 pad = ((TVirtualPad*)l->At(0)); pad->cd();
1231 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
2589cf64 1232 selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1233 if(!GetGraphArray(xy, kTrackIn, 0, 1, n, selection)) break;
2589cf64 1234 xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 2.5;
1ee39b3a 1235 pad=((TVirtualPad*)l->At(1)); pad->cd();
1236 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
81979445 1237 if(!GetGraphArray(xy, kTrackIn, 1, 1)) break;
1ee39b3a 1238 return kTRUE;
6558fd69 1239 case 9: // kTrackIn z
1ee39b3a 1240 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1241 xy[0] = -1.; xy[1] = -1000.; xy[2] = 1.; xy[3] = 4000.;
1242 pad = ((TVirtualPad*)l->At(0)); pad->cd();
1243 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
dfd7d48b 1244 selection[0]=1;
1245 if(!GetGraphArray(xy, kTrackIn, 2, 1, 1, selection)) break;
1ee39b3a 1246 xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
1247 pad = ((TVirtualPad*)l->At(1)); pad->cd();
1248 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
81979445 1249 if(!GetGraphArray(xy, kTrackIn, 3, 1)) break;
1ee39b3a 1250 return kTRUE;
6558fd69 1251 case 10: // kTrackIn phi
1ee39b3a 1252 xy[0] = -.3; xy[1] = -5.; xy[2] = .3; xy[3] = 50.;
a310e49b 1253 if(GetGraph(&xy[0], kTrackIn, 4)) return kTRUE;
1254 break;
6558fd69 1255 case 11: // kTrackOut y
1256 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1257 xy[0] = -.3; xy[1] = -50.; xy[2] = .3; xy[3] = 150.;
1258 pad = ((TVirtualPad*)l->At(0)); pad->cd();
1259 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
2589cf64 1260 selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1261 if(!GetGraphArray(xy, kTrackOut, 0, 1, n, selection)) break;
1262 pad=((TVirtualPad*)l->At(1)); pad->cd();
1263 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
2589cf64 1264 selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1265 if(!GetGraphArray(xy, kTrackOut, 0, 1, n, selection)) break;
1266 return kTRUE;
1267 case 12: // kTrackOut y
a310e49b 1268 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
6558fd69 1269 xy[0] = -.3; xy[1] = -50.; xy[2] = .3; xy[3] = 150.;
a310e49b 1270 pad = ((TVirtualPad*)l->At(0)); pad->cd();
1271 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
2589cf64 1272 selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1273 if(!GetGraphArray(xy, kTrackOut, 0, 1, n, selection)) break;
2589cf64 1274 xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 2.5;
a310e49b 1275 pad=((TVirtualPad*)l->At(1)); pad->cd();
1276 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
81979445 1277 if(!GetGraphArray(xy, kTrackOut, 1, 1)) break;
a310e49b 1278 return kTRUE;
6558fd69 1279 case 13: // kTrackOut z
a310e49b 1280 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1281 xy[0] = -1.; xy[1] = -1000.; xy[2] = 1.; xy[3] = 4000.;
1282 pad = ((TVirtualPad*)l->At(0)); pad->cd();
1283 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
81979445 1284 if(!GetGraphArray(xy, kTrackOut, 2, 1)) break;
a310e49b 1285 xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
1286 pad = ((TVirtualPad*)l->At(1)); pad->cd();
1287 pad->SetMargin(0.1, 0.1, 0.1, 0.01);
81979445 1288 if(!GetGraphArray(xy, kTrackOut, 3, 1)) break;
a310e49b 1289 return kTRUE;
6558fd69 1290 case 14: // kTrackOut phi
a310e49b 1291 xy[0] = -.3; xy[1] = -5.; xy[2] = .3; xy[3] = 50.;
1292 if(GetGraph(&xy[0], kTrackOut, 4)) return kTRUE;
1ee39b3a 1293 break;
6558fd69 1294 case 15: // kMCcluster
92c40e64 1295 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1296 xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3]=650.;
1297 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1298 selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1299 if(!GetGraphArray(xy, kMCcluster, 0, 1, n, selection)) break;
1300 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1301 selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1302 if(!GetGraphArray(xy, kMCcluster, 0, 1, n, selection)) break;
1303 return kTRUE;
1304 case 16: // kMCcluster
1305 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1306 xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3]=650.;
1307 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1308 selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
6558fd69 1309 if(!GetGraphArray(xy, kMCcluster, 0, 1, n, selection)) break;
92c40e64 1310 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1311 xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
81979445 1312 if(!GetGraphArray(xy, kMCcluster, 1, 1)) break;
92c40e64 1313 return kTRUE;
6558fd69 1314 case 17: //kMCtracklet [y]
1ee39b3a 1315 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
92d6d80c 1316 xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =500.;
1ee39b3a 1317 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1318 selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
d25116b6 1319 if(!GetGraphArray(xy, kMCtracklet, 0, 1, n, selection)) break;
1320 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1321 selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
d25116b6 1322 if(!GetGraphArray(xy, kMCtracklet, 0, 1, n, selection)) break;
1323 return kTRUE;
1324 case 18: //kMCtracklet [y]
1325 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1326 xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =500.;
1327 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1328 selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
d25116b6 1329 if(!GetGraphArray(xy, kMCtracklet, 0, 1, n, selection)) break;
1ee39b3a 1330 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1331 xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
81979445 1332 if(!GetGraphArray(xy, kMCtracklet, 1, 1)) break;
1ee39b3a 1333 return kTRUE;
d25116b6 1334 case 19: //kMCtracklet [z]
92c40e64 1335 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1336 xy[0]=-1.; xy[1]=-100.; xy[2]=1.; xy[3] =2500.;
1337 ((TVirtualPad*)l->At(0))->cd();
81979445 1338 if(!GetGraphArray(xy, kMCtracklet, 2)) break;
92c40e64 1339 xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
1340 ((TVirtualPad*)l->At(1))->cd();
81979445 1341 if(!GetGraphArray(xy, kMCtracklet, 3)) break;
92c40e64 1342 return kTRUE;
d25116b6 1343 case 20: //kMCtracklet [phi]
92c40e64 1344 xy[0]=-.3; xy[1]=-3.; xy[2]=.3; xy[3] =25.;
a310e49b 1345 if(!GetGraph(&xy[0], kMCtracklet, 4)) break;
1ee39b3a 1346 return kTRUE;
d25116b6 1347 case 21: //kMCtrack [y] ly [0]
1348 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1349 xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
1350 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1351 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*0.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1352 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer1")) break;
1353 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1354 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*0.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1355 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer1")) break;
1356 return kTRUE;
1357 case 22: //kMCtrack [y] ly [1]
1ee39b3a 1358 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1295b37f 1359 xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
1ee39b3a 1360 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1361 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*1.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1362 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer2")) break;
1363 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1364 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*1.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1365 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer2")) break;
1366 return kTRUE;
1367 case 23: //kMCtrack [y] ly [2]
1368 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1369 xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
1370 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1371 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*2.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1372 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer3")) break;
1373 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1374 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*2.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1375 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer3")) break;
1376 return kTRUE;
1377 case 24: //kMCtrack [y] ly [3]
1378 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1379 xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
1380 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1381 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*3.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1382 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer4")) break;
1383 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1384 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*3.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1385 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer4")) break;
1386 return kTRUE;
1387 case 25: //kMCtrack [y] ly [4]
1388 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1389 xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
1390 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1391 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*4.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1392 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer5")) break;
1393 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1394 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*4.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1395 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer5")) break;
1396 return kTRUE;
1397 case 26: //kMCtrack [y] ly [5]
1398 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1399 xy[0]=-.2; xy[1]=-50.; xy[2]=.2; xy[3] =400.;
1400 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1401 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*5.); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1402 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer6")) break;
1403 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1404 selStart=Int_t(fgkNresYsegm[fSegmentLevel]*5.5); for(n=0; n<fgkNresYsegm[fSegmentLevel]/2; n++) selection[n]=selStart+n;
d25116b6 1405 if(!GetGraphArray(xy, kMCtrack, 0, 1, n, selection, "Layer6")) break;
1406 return kTRUE;
1407 case 27: //kMCtrack [y pulls]
1408 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
2589cf64 1409 xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 5.5;
d25116b6 1410 ((TVirtualPad*)l->At(0))->cd();
81979445 1411 selStart=0; for(n=0; n<6; n++) selection[n]=selStart+n;
d25116b6 1412 if(!GetGraphArray(xy, kMCtrack, 1, 1, n, selection)) break;
1ee39b3a 1413 ((TVirtualPad*)l->At(1))->cd();
81979445 1414 selStart=6; for(n=0; n<6; n++) selection[n]=selStart+n;
d25116b6 1415 if(!GetGraphArray(xy, kMCtrack, 1, 1, n, selection)) break;
1ee39b3a 1416 return kTRUE;
d25116b6 1417 case 28: //kMCtrack [z]
1ee39b3a 1418 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
a310e49b 1419 xy[0]=-1.; xy[1]=-1500.; xy[2]=1.; xy[3] =6000.;
1ee39b3a 1420 ((TVirtualPad*)l->At(0))->cd();
92d6d80c 1421 if(!GetGraphArray(xy, kMCtrack, 2)) break;
a310e49b 1422 xy[0] = -1.; xy[1] = -1.5; xy[2] = 1.; xy[3] = 5.;
1ee39b3a 1423 ((TVirtualPad*)l->At(1))->cd();
92d6d80c 1424 if(!GetGraphArray(xy, kMCtrack, 3)) break;
1ee39b3a 1425 return kTRUE;
d25116b6 1426 case 29: //kMCtrack [phi/snp]
1ee39b3a 1427 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
92d6d80c 1428 xy[0]=-.2; xy[1]=-0.5; xy[2]=.2; xy[3] =10.;
1ee39b3a 1429 ((TVirtualPad*)l->At(0))->cd();
92d6d80c 1430 if(!GetGraphArray(xy, kMCtrack, 4)) break;
a310e49b 1431 xy[0] = -.2; xy[1] = -1.5; xy[2] = .2; xy[3] = 5.;
1ee39b3a 1432 ((TVirtualPad*)l->At(1))->cd();
92d6d80c 1433 if(!GetGraphArray(xy, kMCtrack, 5)) break;
1ee39b3a 1434 return kTRUE;
d25116b6 1435 case 30: //kMCtrack [theta/tgl]
1ee39b3a 1436 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1437 xy[0]=-1.; xy[1]=-0.5; xy[2]=1.; xy[3] =5.;
1438 ((TVirtualPad*)l->At(0))->cd();
92d6d80c 1439 if(!GetGraphArray(xy, kMCtrack, 6)) break;
1ee39b3a 1440 xy[0] = -.2; xy[1] = -0.5; xy[2] = .2; xy[3] = 2.5;
1441 ((TVirtualPad*)l->At(1))->cd();
92d6d80c 1442 if(!GetGraphArray(xy, kMCtrack, 7)) break;
1ee39b3a 1443 return kTRUE;
d25116b6 1444 case 31: //kMCtrack [pt]
92c40e64 1445 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1446 pad = (TVirtualPad*)l->At(0); pad->cd();
1295b37f 1447 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1448 // pi selection
1449 n=0;
1450 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1451 selection[n++] = il*11 + 2; // pi-
1452 selection[n++] = il*11 + 8; // pi+
1453 }
a310e49b 1454 xy[0] = 0.2; xy[1] = -.7; xy[2] = 7.; xy[3] = 4.;
1455 //xy[0] = 0.2; xy[1] = -1.; xy[2] = 7.; xy[3] = 10.; // SA
92d6d80c 1456 if(!GetGraphArray(xy, kMCtrack, 8, kTRUE, n, selection, "#pi#pm")) break;
31c8fa8a 1457 pad->Modified(); pad->Update(); pad->SetLogx();
92c40e64 1458 pad = (TVirtualPad*)l->At(1); pad->cd();
1295b37f 1459 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1460 // mu selection
1461 n=0;
1462 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1463 selection[n++] = il*11 + 3; // mu-
1464 selection[n++] = il*11 + 7; // mu+
1465 }
1466 if(!GetGraphArray(xy, kMCtrack, 8, kTRUE, n, selection, "#mu#pm")) break;
31c8fa8a 1467 pad->Modified(); pad->Update(); pad->SetLogx();
92c40e64 1468 return kTRUE;
d25116b6 1469 case 32: //kMCtrack [pt]
92c40e64 1470 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1471 pad = (TVirtualPad*)l->At(0); pad->cd();
1295b37f 1472 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1473 // p selection
1474 n=0;
1475 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1476 selection[n++] = il*11 + 0; // p bar
1477 selection[n++] = il*11 + 10; // p
1478 }
a310e49b 1479 xy[0] = 0.2; xy[1] = -.7; xy[2] = 7.; xy[3] = 8.;
1480 //xy[0] = 0.2; xy[1] = -1.; xy[2] = 7.; xy[3] = 10.; // SA
92d6d80c 1481 if(!GetGraphArray(xy, kMCtrack, 8, kTRUE, n, selection, "p&p bar")) break;
1482 pad->Modified(); pad->Update(); pad->SetLogx();
92c40e64 1483 pad = (TVirtualPad*)l->At(1); pad->cd();
1295b37f 1484 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1485 // e selection
1486 n=0;
1487 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1488 selection[n++] = il*11 + 4; // e-
1489 selection[n++] = il*11 + 6; // e+
1490 }
a310e49b 1491 xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 12.;
1492 //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 14.; // SA
92d6d80c 1493 if(!GetGraphArray(xy, kMCtrack, 8, kTRUE, n, selection, "e#pm")) break;
1494 pad->Modified(); pad->Update(); pad->SetLogx();
92c40e64 1495 return kTRUE;
d25116b6 1496 case 33: //kMCtrack [1/pt] pulls
a310e49b 1497 xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5;
1498 //xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 4.5; // SA
92c40e64 1499 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1500 pad = (TVirtualPad*)l->At(0); pad->cd();
1295b37f 1501 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1502 // pi selection
1503 n=0;
1504 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1505 selection[n++] = il*11 + 2; // pi-
1506 selection[n++] = il*11 + 8; // pi+
1507 }
1508 if(!GetGraphArray(xy, kMCtrack, 9, kTRUE, n, selection, "#pi#pm")) break;
92c40e64 1509 pad = (TVirtualPad*)l->At(1); pad->cd();
1295b37f 1510 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1511 // mu selection
1512 n=0;
1513 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1514 selection[n++] = il*11 + 3; // mu-
1515 selection[n++] = il*11 + 7; // mu+
1516 }
1517 if(!GetGraphArray(xy, kMCtrack, 9, kTRUE, n, selection, "#mu#pm")) break;
92c40e64 1518 return kTRUE;
d25116b6 1519 case 34: //kMCtrack [1/pt] pulls
92c40e64 1520 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1521 pad = (TVirtualPad*)l->At(0); pad->cd();
1295b37f 1522 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1523 // p selection
1524 n=0;
1525 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1526 selection[n++] = il*11 + 0; // p bar
1527 selection[n++] = il*11 + 10; // p
1528 }
a310e49b 1529 xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5;
1530 //xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 6.; // SA
92d6d80c 1531 if(!GetGraphArray(xy, kMCtrack, 9, kTRUE, n, selection, "p & p bar")) break;
92c40e64 1532 pad = (TVirtualPad*)l->At(1); pad->cd();
1295b37f 1533 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1534 // e selection
1535 n=0;
1536 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1537 selection[n++] = il*11 + 4; // e-
1538 selection[n++] = il*11 + 6; // e+
1539 }
a310e49b 1540 xy[0] = 0.; xy[1] = -2.; xy[2] = 2.; xy[3] = 4.5;
92d6d80c 1541 if(!GetGraphArray(xy, kMCtrack, 9, kTRUE, n, selection, "e#pm")) break;
1ee39b3a 1542 return kTRUE;
d25116b6 1543 case 35: //kMCtrack [p]
a310e49b 1544 xy[0] = 0.2; xy[1] = -.7; xy[2] = 7.; xy[3] = 4.;
1545 //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 10.;
92c40e64 1546 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1547 pad = (TVirtualPad*)l->At(0); pad->cd();
1295b37f 1548 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1549 // pi selection
1550 n=0;
1551 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1552 selection[n++] = il*11 + 2; // pi-
1553 selection[n++] = il*11 + 8; // pi+
1554 }
1555 if(!GetGraphArray(xy, kMCtrack, 10, kTRUE, n, selection, "#pi#pm")) break;
1556 pad->Modified(); pad->Update(); pad->SetLogx();
92c40e64 1557 pad = (TVirtualPad*)l->At(1); pad->cd();
1295b37f 1558 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1559 // mu selection
1560 n=0;
1561 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1562 selection[n++] = il*11 + 3; // mu-
1563 selection[n++] = il*11 + 7; // mu+
1564 }
1565 if(!GetGraphArray(xy, kMCtrack, 10, kTRUE, n, selection, "#mu#pm")) break;
1566 pad->Modified(); pad->Update(); pad->SetLogx();
1ee39b3a 1567 return kTRUE;
d25116b6 1568 case 36: //kMCtrack [p]
92c40e64 1569 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1570 pad = (TVirtualPad*)l->At(0); pad->cd();
1295b37f 1571 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1572 // p selection
1573 n=0;
1574 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1575 selection[n++] = il*11 + 0; // p bar
1576 selection[n++] = il*11 + 10; // p
1577 }
a310e49b 1578 xy[0] = 0.2; xy[1] = -.7; xy[2] = 7.; xy[3] = 8.;
1579 //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 12.; // SA
92d6d80c 1580 if(!GetGraphArray(xy, kMCtrack, 10, kTRUE, n, selection, "p & p bar")) break;
1581 pad->Modified(); pad->Update(); pad->SetLogx();
92c40e64 1582 pad = (TVirtualPad*)l->At(1); pad->cd();
1295b37f 1583 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
92d6d80c 1584 // e selection
1585 n=0;
1586 for(Int_t il(ly0); il<AliTRDgeometry::kNlayer; il+=dly){
1587 selection[n++] = il*11 + 4; // e-
1588 selection[n++] = il*11 + 6; // e+
1589 }
a310e49b 1590 xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 12.;
1591 //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 14.; // SA
92d6d80c 1592 if(!GetGraphArray(xy, kMCtrack, 10, kTRUE, n, selection, "e#pm")) break;
1593 pad->Modified(); pad->Update(); pad->SetLogx();
1ee39b3a 1594 return kTRUE;
d25116b6 1595 case 37: // kMCtrackIn [y]
1ee39b3a 1596 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
a310e49b 1597 xy[0]=-.25; xy[1]=-1000.; xy[2]=.25; xy[3] =3000.;
1ee39b3a 1598 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1599 selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
a310e49b 1600 if(!GetGraphArray(xy, kMCtrackIn, 0, 1, n, selection)) break;
1601 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1602 selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
a310e49b 1603 if(!GetGraphArray(&xy[0], kMCtrackIn, 0, 1, n, selection)) break;
1604 return kTRUE;
d25116b6 1605 case 38: // kMCtrackIn [y]
a310e49b 1606 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1607 xy[0]=-.25; xy[1]=-1000.; xy[2]=.25; xy[3] =3000.;
1608 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1609 selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
a310e49b 1610 if(!GetGraphArray(xy, kMCtrackIn, 0, 1, n, selection)) break;
2589cf64 1611 xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 2.5;
1ee39b3a 1612 ((TVirtualPad*)l->At(1))->cd();
81979445 1613 if(!GetGraphArray(xy, kMCtrackIn, 1, 1)) break;
1ee39b3a 1614 return kTRUE;
d25116b6 1615 case 39: // kMCtrackIn [z]
1ee39b3a 1616 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1617 xy[0]=-1.; xy[1]=-500.; xy[2]=1.; xy[3] =800.;
1618 ((TVirtualPad*)l->At(0))->cd();
81979445 1619 if(!GetGraphArray(xy, kMCtrackIn, 2, 1)) break;
1ee39b3a 1620 xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
1621 ((TVirtualPad*)l->At(1))->cd();
81979445 1622 if(!GetGraphArray(xy, kMCtrackIn, 3, 1)) break;
1ee39b3a 1623 return kTRUE;
d25116b6 1624 case 40: // kMCtrackIn [phi|snp]
1ee39b3a 1625 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1626 xy[0]=-.25; xy[1]=-0.5; xy[2]=.25; xy[3] =2.5;
1627 ((TVirtualPad*)l->At(0))->cd();
a310e49b 1628 if(!GetGraph(&xy[0], kMCtrackIn, 4)) break;
1ee39b3a 1629 xy[0] = -.25; xy[1] = -0.5; xy[2] = .25; xy[3] = 1.5;
1630 ((TVirtualPad*)l->At(1))->cd();
a310e49b 1631 if(!GetGraph(&xy[0], kMCtrackIn, 5)) break;
1ee39b3a 1632 return kTRUE;
d25116b6 1633 case 41: // kMCtrackIn [theta|tgl]
1ee39b3a 1634 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1635 xy[0]=-1.; xy[1]=-1.; xy[2]=1.; xy[3] =4.;
1636 ((TVirtualPad*)l->At(0))->cd();
a310e49b 1637 if(!GetGraph(&xy[0], kMCtrackIn, 6)) break;
1ee39b3a 1638 xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 1.5;
1639 ((TVirtualPad*)l->At(1))->cd();
a310e49b 1640 if(!GetGraph(&xy[0], kMCtrackIn, 7)) break;
1ee39b3a 1641 return kTRUE;
d25116b6 1642 case 42: // kMCtrackIn [pt]
92d6d80c 1643 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
a310e49b 1644 xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.;
1645 //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 10.; // SA
92d6d80c 1646 pad=(TVirtualPad*)l->At(0); pad->cd(); pad->SetLogx();
1647 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
1648 n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
1649 if(!GetGraphArray(xy, kMCtrackIn, 8, 1, n, selection)) break;
1650 pad = (TVirtualPad*)l->At(1); pad->cd(); pad->SetLogx();
1651 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
1652 n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
1653 if(!GetGraphArray(xy, kMCtrackIn, 8, 1, n, selection)) break;
1654 return kTRUE;
d25116b6 1655 case 43: //kMCtrackIn [1/pt] pulls
92d6d80c 1656 xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5;
1657 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1658 pad = (TVirtualPad*)l->At(0); pad->cd();
1659 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
1660 n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
1661 if(!GetGraphArray(xy, kMCtrackIn, 9, 1, n, selection)) break;
1662 pad = (TVirtualPad*)l->At(1); pad->cd();
1663 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
1664 n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
1665 if(!GetGraphArray(xy, kMCtrackIn, 9, 1, n, selection)) break;
1666 return kTRUE;
d25116b6 1667 case 44: // kMCtrackIn [p]
a310e49b 1668 xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.;
1669 //xy[0] = 0.2; xy[1] = -1.5; xy[2] = 7.; xy[3] = 10.;
92d6d80c 1670 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1671 pad = ((TVirtualPad*)l->At(0));pad->cd();pad->SetLogx();
1672 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
1673 n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
1674 if(!GetGraphArray(xy, kMCtrackIn, 10, 1, n, selection)) break;
1675 pad = ((TVirtualPad*)l->At(1)); pad->cd();pad->SetLogx();
1676 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
1677 n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
1678 if(!GetGraphArray(xy, kMCtrackIn, 10, 1, n, selection)) break;
1679 return kTRUE;
d25116b6 1680 case 45: // kMCtrackOut [y]
92d6d80c 1681 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
a310e49b 1682 xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =400.;
92d6d80c 1683 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1684 selStart=0; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
a310e49b 1685 if(!GetGraphArray(xy, kMCtrackOut, 0, 1, n, selection)) break;
1686 ((TVirtualPad*)l->At(1))->cd();
2589cf64 1687 selStart=fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
a310e49b 1688 if(!GetGraphArray(&xy[0], kMCtrackOut, 0, 1, n, selection)) break;
1689 return kTRUE;
d25116b6 1690 case 46: // kMCtrackOut [y]
a310e49b 1691 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1692 xy[0]=-.3; xy[1]=-50.; xy[2]=.3; xy[3] =400.;
1693 ((TVirtualPad*)l->At(0))->cd();
2589cf64 1694 selStart=2*fgkNresYsegm[fSegmentLevel]/3; for(n=0; n<fgkNresYsegm[fSegmentLevel]/3; n++) selection[n]=selStart+n;
a310e49b 1695 if(!GetGraphArray(xy, kMCtrackOut, 0, 1, n, selection)) break;
2589cf64 1696 xy[0] = -.5; xy[1] = -0.5; xy[2] = fgkNresYsegm[fSegmentLevel]-.5; xy[3] = 2.5;
92d6d80c 1697 ((TVirtualPad*)l->At(1))->cd();
81979445 1698 if(!GetGraphArray(xy, kMCtrackOut, 1, 1)) break;
92d6d80c 1699 return kTRUE;
d25116b6 1700 case 47: // kMCtrackOut [z]
92d6d80c 1701 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
a310e49b 1702 xy[0]=-1.; xy[1]=-500.; xy[2]=1.; xy[3] =1500.;
92d6d80c 1703 ((TVirtualPad*)l->At(0))->cd();
81979445 1704 if(!GetGraphArray(xy, kMCtrackOut, 2, 1)) break;
92d6d80c 1705 xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 2.5;
1706 ((TVirtualPad*)l->At(1))->cd();
81979445 1707 if(!GetGraphArray(xy, kMCtrackOut, 3, 1)) break;
92d6d80c 1708 return kTRUE;
d25116b6 1709 case 48: // kMCtrackOut [phi|snp]
92d6d80c 1710 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1711 xy[0]=-.25; xy[1]=-0.5; xy[2]=.25; xy[3] =2.5;
1712 ((TVirtualPad*)l->At(0))->cd();
a310e49b 1713 if(!GetGraph(&xy[0], kMCtrackOut, 4)) break;
92d6d80c 1714 xy[0] = -.25; xy[1] = -0.5; xy[2] = .25; xy[3] = 1.5;
1715 ((TVirtualPad*)l->At(1))->cd();
a310e49b 1716 if(!GetGraph(&xy[0], kMCtrackOut, 5)) break;
92d6d80c 1717 return kTRUE;
d25116b6 1718 case 49: // kMCtrackOut [theta|tgl]
92d6d80c 1719 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1720 xy[0]=-1.; xy[1]=-1.; xy[2]=1.; xy[3] =4.;
1721 ((TVirtualPad*)l->At(0))->cd();
a310e49b 1722 if(!GetGraph(&xy[0], kMCtrackOut, 6)) break;
1723 xy[0] = -1.; xy[1] = -0.5; xy[2] = 1.; xy[3] = 15.;
92d6d80c 1724 ((TVirtualPad*)l->At(1))->cd();
a310e49b 1725 if(!GetGraph(&xy[0], kMCtrackOut, 7)) break;
92d6d80c 1726 return kTRUE;
d25116b6 1727 case 50: // kMCtrackOut [pt]
1ee39b3a 1728 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
db99a57a 1729 xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.;
a310e49b 1730 pad=(TVirtualPad*)l->At(0); pad->cd(); pad->SetLogx();
1295b37f 1731 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
a310e49b 1732 n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
1733 if(!GetGraphArray(xy, kMCtrackOut, 8, 1, n, selection)) break;
1734 pad = (TVirtualPad*)l->At(1); pad->cd();pad->SetLogx();
1295b37f 1735 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
a310e49b 1736 n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
1737 if(!GetGraphArray(xy, kMCtrackOut, 8, 1, n, selection)) break;
1ee39b3a 1738 return kTRUE;
d25116b6 1739 case 51: //kMCtrackOut [1/pt] pulls
92d6d80c 1740 xy[0] = 0.; xy[1] = -1.; xy[2] = 2.; xy[3] = 3.5;
1741 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
1742 pad = (TVirtualPad*)l->At(0); pad->cd();
1743 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
a310e49b 1744 n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
1745 if(!GetGraphArray(xy, kMCtrackOut, 9, 1, n, selection)) break;
92d6d80c 1746 pad = (TVirtualPad*)l->At(1); pad->cd();
1747 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
a310e49b 1748 n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
1749 if(!GetGraphArray(xy, kMCtrackOut, 9, 1, n, selection)) break;
92d6d80c 1750 return kTRUE;
d25116b6 1751 case 52: // kMCtrackOut [p]
1ee39b3a 1752 gPad->Divide(2, 1, 1.e-5, 1.e-5); l=gPad->GetListOfPrimitives();
db99a57a 1753 xy[0] = 0.2; xy[1] = -.8; xy[2] = 7.; xy[3] = 6.;
a310e49b 1754 pad = ((TVirtualPad*)l->At(0));pad->cd();pad->SetLogx();
1295b37f 1755 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
a310e49b 1756 n=0; selection[n++]=2; selection[n++]=3; selection[n++]=7; selection[n++]=8;
1757 if(!GetGraphArray(xy, kMCtrackOut, 10, 1, n, selection)) break;
1758 pad = ((TVirtualPad*)l->At(1)); pad->cd();pad->SetLogx();
1295b37f 1759 pad->SetMargin(0.125, 0.015, 0.1, 0.015);
a310e49b 1760 n=0; selection[n++]=0; selection[n++]=4; selection[n++]=6; selection[n++]=10;
1761 if(!GetGraphArray(xy, kMCtrackOut, 10, 1, n, selection)) break;
1ee39b3a 1762 return kTRUE;
1763 }
1764 AliWarning(Form("Reference plot [%d] missing result", ifig));
1765 return kFALSE;
1766}
1767
00a3f7a4 1768//________________________________________________________
1769void AliTRDresolution::MakeSummary()
1770{
1771// Build summary plots
1772
1773 if(!fGraphS || !fGraphM){
1774 AliError("Missing results");
1775 return;
1776 }
1777 Float_t xy[4] = {0., 0., 0., 0.};
068e2c00 1778 Float_t range[2];
00a3f7a4 1779 TH2 *h2 = new TH2I("h2SF", "", 20, -.2, .2, fgkNresYsegm[fSegmentLevel], -0.5, fgkNresYsegm[fSegmentLevel]-0.5);
1780 h2->GetXaxis()->CenterTitle();
1781 h2->GetYaxis()->CenterTitle();
068e2c00 1782 h2->GetZaxis()->CenterTitle();h2->GetZaxis()->SetTitleOffset(1.4);
1783
1784 Int_t ih2(0), iSumPlot(0);
1785 TCanvas *cOut = new TCanvas(Form("TRDsummary%s_%d", GetName(), iSumPlot++), "Cluster & Tracklet Resolution", 1024, 768);
1786 cOut->Divide(3,2, 2.e-3, 2.e-3, kYellow-7);
1787 TVirtualPad *p(NULL);
1788
1789 p=cOut->cd(1);
1790 p->SetRightMargin(0.16);p->SetTopMargin(0.06);
1791 h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
1792 h2->SetTitle(Form("Cluster-Track R-Phi Residuals;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel]));
1793 MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphS->At(kCluster))->At(0), h2);
1794 GetRange(h2, 1, range);
1795 h2->GetZaxis()->SetRangeUser(range[0], range[1]);
1796 h2->Draw("colz");
1797 h2->SetContour(7);
1798
1799 p=cOut->cd(2);
1800 p->SetRightMargin(0.16);p->SetTopMargin(0.06);
1801 h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
1802 h2->SetTitle(Form("Cluster-Track R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel]));
1803 MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphM->At(kCluster))->At(0), h2);
1804 GetRange(h2, 0, range);
1805 h2->GetZaxis()->SetRangeUser(range[0], range[1]);
1806 h2->Draw("colz");
1807 h2->SetContour(7);
1808
1809 p=cOut->cd(3);
1810 p->SetRightMargin(0.06);p->SetTopMargin(0.06);
1811 xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
1812 GetGraphArray(xy, kCluster, 1, 1);
1813
1814 p=cOut->cd(4);
1815 p->SetRightMargin(0.16);p->SetTopMargin(0.06);
1816 h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
1817 h2->SetTitle(Form("Tracklet-Track R-Phi Residuals;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel]));
1818 MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphS->At(kTrack))->At(0), h2);
1819 GetRange(h2, 1, range);
1820 h2->GetZaxis()->SetRangeUser(range[0], range[1]);
1821 h2->Draw("colz");
1822 h2->SetContour(7);
1823
1824 p=cOut->cd(5);
1825 p->SetRightMargin(0.16);p->SetTopMargin(0.06);
1826 h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
1827 h2->SetTitle(Form("Tracklet-Track R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel]));
1828 MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphM->At(kTrack))->At(0), h2);
1829 GetRange(h2, 0, range);
1830 h2->GetZaxis()->SetRangeUser(range[0], range[1]);
1831 h2->Draw("colz");
1832 h2->SetContour(7);
1833
1834 p=cOut->cd(6);
1835 p->SetRightMargin(0.06);p->SetTopMargin(0.06);
1836 xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
1837 GetGraphArray(xy, kTrack, 1, 1);
1838
1839 cOut->SaveAs(Form("%s.gif", cOut->GetName()));
00a3f7a4 1840
068e2c00 1841 if(!HasMCdata()){
1842 delete cOut;
1843 return;
1844 }
1845 cOut->Clear(); cOut->SetName(Form("TRDsummary%s_%d", GetName(), iSumPlot++));
1846 cOut->Divide(3, 2, 2.e-3, 2.e-3, kBlue-10);
1847
1848 p=cOut->cd(1);
1849 p->SetRightMargin(0.16);p->SetTopMargin(0.06);
1850 h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
1851 h2->SetTitle(Form("Cluster-MC R-Phi Resolution;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel]));
1852 MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphS->At(kMCcluster))->At(0), h2);
1853 GetRange(h2, 1, range);
1854 h2->GetZaxis()->SetRangeUser(range[0], range[1]);
1855 h2->Draw("colz");
1856 h2->SetContour(7);
1857
1858 p=cOut->cd(2);
1859 p->SetRightMargin(0.16);p->SetTopMargin(0.06);
1860 h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
1861 h2->SetContour(7);
1862 h2->SetTitle(Form("Cluster-MC R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel]));
1863 MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphM->At(kMCcluster))->At(0), h2);
1864 GetRange(h2, 0, range);
1865 h2->GetZaxis()->SetRangeUser(range[0], range[1]);
1866 h2->Draw("colz");
1867 h2->SetContour(7);
1868
1869 p=cOut->cd(3);
1870 p->SetRightMargin(0.06);p->SetTopMargin(0.06);
00a3f7a4 1871 xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
1872 GetGraphArray(xy, kMCcluster, 1, 1);
1873
068e2c00 1874 p=cOut->cd(4);
1875 p->SetRightMargin(0.16);p->SetTopMargin(0.06);
1876 h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
1877 h2->SetContour(7);
1878 h2->SetTitle(Form("Tracklet-MC R-Phi Resolution;tg(#phi);%s;Sigma [#mum]", fgkResYsegmName[fSegmentLevel]));
1879 MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphS->At(kMCtracklet))->At(0), h2);
1880 GetRange(h2, 1, range);
1881 h2->GetZaxis()->SetRangeUser(range[0], range[1]);
1882 h2->Draw("colz");
1883 h2->SetContour(7);
1884
1885 p=cOut->cd(5);
1886 p->SetRightMargin(0.16);p->SetTopMargin(0.06);
1887 h2=(TH2I*)h2->Clone(Form("h2SF_%d", ih2++));
1888 h2->SetContour(7);
1889 h2->SetTitle(Form("Tracklet-MC R-Phi Systematics;tg(#phi);%s;Mean [#mum]", fgkResYsegmName[fSegmentLevel]));
1890 MakeSummaryPlot((TObjArray*) ((TObjArray*)fGraphM->At(kMCtracklet))->At(0), h2);
1891 GetRange(h2, 0, range);
1892 h2->GetZaxis()->SetRangeUser(range[0], range[1]);
1893 h2->Draw("colz");
1894 h2->SetContour(7);
1895
1896 p=cOut->cd(6);
1897 p->SetRightMargin(0.06);p->SetTopMargin(0.06);
1898 xy[0]=-.5; xy[1]=-0.5; xy[2]=fgkNresYsegm[fSegmentLevel]-.5; xy[3]=2.5;
1899 GetGraphArray(xy, kMCtracklet, 1, 1);
1900
1901 cOut->SaveAs(Form("%s.gif", cOut->GetName()));
1902 delete cOut;
1903}
1904
1905//________________________________________________________
1906void AliTRDresolution::GetRange(TH2 *h2, Char_t mod, Float_t *range)
1907{
1908// Returns the range of the bulk of data in histogram h2. Removes outliers.
1909// The "range" vector should be initialized with 2 elements
1910// Option "mod" can be any of
1911// - 0 : gaussian like distribution
1912// - 1 : tailed distribution
1913
1914 Int_t nx(h2->GetNbinsX())
1915 , ny(h2->GetNbinsY())
1916 , n(nx*ny);
1917 Double_t *data=new Double_t[n];
1918 for(Int_t ix(1), in(0); ix<=nx; ix++){
1919 for(Int_t iy(1); iy<=ny; iy++)
1920 data[in++] = h2->GetBinContent(ix, iy);
00a3f7a4 1921 }
068e2c00 1922 Double_t mean, sigm;
1923 AliMathBase::EvaluateUni(n, data, mean, sigm, Int_t(n*.8));
1924
1925 range[0]=mean-3.*sigm; range[1]=mean+3.*sigm;
1926 if(mod==1) range[0]=TMath::Max(Float_t(1.e-3), range[0]);
1927 AliDebug(2, Form("h[%s] range0[%f %f]", h2->GetName(), range[0], range[1]));
1928 TH1S h1("h1SF0", "", 100, range[0], range[1]);
1929 h1.FillN(n,data,0);
1930 delete [] data;
1931
1932 switch(mod){
1933 case 0:// gaussian distribution
1934 {
1935 TF1 fg("fg", "gaus", mean-3.*sigm, mean+3.*sigm);
1936 h1.Fit(&fg, "QN");
1937 mean=fg.GetParameter(1); sigm=fg.GetParameter(2);
1938 range[0] = mean-2.5*sigm;range[1] = mean+2.5*sigm;
1939 AliDebug(2, Form(" rangeG[%f %f]", range[0], range[1]));
1940 break;
00a3f7a4 1941 }
068e2c00 1942 case 1:// tailed distribution
1943 {
1944 Int_t bmax(h1.GetMaximumBin());
1945 Int_t jBinMin(1), jBinMax(100);
1946 for(Int_t ibin(bmax); ibin--;){
1947 if(h1.GetBinContent(ibin)==0){
1948 jBinMin=ibin; break;
1949 }
1950 }
1951 for(Int_t ibin(bmax); ibin++;){
1952 if(h1.GetBinContent(ibin)==0){
1953 jBinMax=ibin; break;
1954 }
00a3f7a4 1955 }
068e2c00 1956 range[0]=h1.GetBinCenter(jBinMin); range[1]=h1.GetBinCenter(jBinMax);
1957 AliDebug(2, Form(" rangeT[%f %f]", range[0], range[1]));
1958 break;
1959 }
00a3f7a4 1960 }
00a3f7a4 1961
068e2c00 1962 return;
1963}
1964
1965//________________________________________________________
1966void AliTRDresolution::MakeSummaryPlot(TObjArray *a, TH2 *h2)
1967{
1968 h2->Reset();
1969 Double_t x,y;
1970 TGraphErrors *g(NULL); TAxis *ax(h2->GetXaxis());
00a3f7a4 1971 for(Int_t iseg(0); iseg<fgkNresYsegm[fSegmentLevel]; iseg++){
1972 g=(TGraphErrors*)a->At(iseg);
1973 for(Int_t in(0); in<g->GetN(); in++){
1974 g->GetPoint(in, x, y);
1975 h2->SetBinContent(ax->FindBin(x), iseg+1, y);
00a3f7a4 1976 }
1977 }
00a3f7a4 1978}
1979
068e2c00 1980
00a3f7a4 1981//________________________________________________________
92c40e64 1982Char_t const *fgParticle[11]={
1983 " p bar", " K -", " #pi -", " #mu -", " e -",
1984 " No PID",
1985 " e +", " #mu +", " #pi +", " K +", " p",
1986};
1295b37f 1987const Color_t fgColorS[11]={
1988kOrange, kOrange-3, kMagenta+1, kViolet, kRed,
31c8fa8a 1989kGray,
1295b37f 1990kRed, kViolet, kMagenta+1, kOrange-3, kOrange
1991};
1992const Color_t fgColorM[11]={
1993kCyan-5, kAzure-4, kBlue-7, kBlue+2, kViolet+10,
1994kBlack,
1995kViolet+10, kBlue+2, kBlue-7, kAzure-4, kCyan-5
31c8fa8a 1996};
1997const Marker_t fgMarker[11]={
199830, 30, 26, 25, 24,
199928,
200020, 21, 22, 29, 29
2001};
1ee39b3a 2002Bool_t AliTRDresolution::PostProcess()
2003{
2004 //fContainer = dynamic_cast<TObjArray*>(GetOutputData(0));
2005 if (!fContainer) {
2006 AliError("ERROR: list not available");
2007 return kFALSE;
2008 }
b91fdd71 2009 TGraph *gm= NULL, *gs= NULL;
1ee39b3a 2010 if(!fGraphS && !fGraphM){
31c8fa8a 2011 TObjArray *aM(NULL), *aS(NULL);
1ee39b3a 2012 Int_t n = fContainer->GetEntriesFast();
2013 fGraphS = new TObjArray(n); fGraphS->SetOwner();
2014 fGraphM = new TObjArray(n); fGraphM->SetOwner();
92d6d80c 2015 for(Int_t ig(0), nc(0); ig<n; ig++){
92c40e64 2016 fGraphM->AddAt(aM = new TObjArray(fgNproj[ig]), ig);
2017 fGraphS->AddAt(aS = new TObjArray(fgNproj[ig]), ig);
92d6d80c 2018
2019 for(Int_t ic=0; ic<fgNproj[ig]; ic++, nc++){
2589cf64 2020 AliDebug(2, Form("building G[%d] P[%d] N[%d]", ig, ic, fNcomp[nc]));
2021 if(fNcomp[nc]>1){
31c8fa8a 2022 TObjArray *agS(NULL), *agM(NULL);
2589cf64 2023 aS->AddAt(agS = new TObjArray(fNcomp[nc]), ic);
2024 aM->AddAt(agM = new TObjArray(fNcomp[nc]), ic);
2025 for(Int_t is=fNcomp[nc]; is--;){
31c8fa8a 2026 agS->AddAt(gs = new TGraphErrors(), is);
92d6d80c 2027 Int_t is0(is%11), il0(is/11);
31c8fa8a 2028 gs->SetMarkerStyle(fgMarker[is0]);
1295b37f 2029 gs->SetMarkerColor(fgColorS[is0]);
2030 gs->SetLineColor(fgColorS[is0]);
92d6d80c 2031 gs->SetLineStyle(il0);gs->SetLineWidth(2);
2589cf64 2032 gs->SetName(Form("s_%d_%02d_%02d", ig, ic, is));
31c8fa8a 2033
2034 agM->AddAt(gm = new TGraphErrors(), is);
2035 gm->SetMarkerStyle(fgMarker[is0]);
1295b37f 2036 gm->SetMarkerColor(fgColorM[is0]);
2037 gm->SetLineColor(fgColorM[is0]);
92d6d80c 2038 gm->SetLineStyle(il0);gm->SetLineWidth(2);
2589cf64 2039 gm->SetName(Form("m_%d_%02d_%02d", ig, ic, is));
a310e49b 2040 // this is important for labels in the legend
2041 if(ic==0) {
2589cf64 2042 gs->SetTitle(Form("%s %02d", fgkResYsegmName[fSegmentLevel], is%fgkNresYsegm[fSegmentLevel]));
2043 gm->SetTitle(Form("%s %02d", fgkResYsegmName[fSegmentLevel], is%fgkNresYsegm[fSegmentLevel]));
81979445 2044 } else if(ic==1) {
2045 gs->SetTitle(Form("%s Ly[%d]", is%2 ?"z":"y", is/2));
2046 gm->SetTitle(Form("%s Ly[%d]", is%2?"z":"y", is/2));
2047 } else if(ic==2||ic==3) {
2048 gs->SetTitle(Form("%s Ly[%d]", is%2 ?"RC":"no RC", is/2));
2049 gm->SetTitle(Form("%s Ly[%d]", is%2?"RC":"no RC", is/2));
a310e49b 2050 } else if(ic<=7) {
d25116b6 2051 gs->SetTitle(Form("Layer[%d]", is%AliTRDgeometry::kNlayer));
2052 gm->SetTitle(Form("Layer[%d]", is%AliTRDgeometry::kNlayer));
a310e49b 2053 } else {
2054 gs->SetTitle(Form("%s @ ly[%d]", fgParticle[is0], il0));
2055 gm->SetTitle(Form("%s @ ly[%d]", fgParticle[is0], il0));
2056 }
92c40e64 2057 }
2058 } else {
2059 aS->AddAt(gs = new TGraphErrors(), ic);
2060 gs->SetMarkerStyle(23);
2061 gs->SetMarkerColor(kRed);
2062 gs->SetLineColor(kRed);
2589cf64 2063 gs->SetNameTitle(Form("s_%d_%02d", ig, ic), "sigma");
92c40e64 2064
2065 aM->AddAt(gm = ig ? (TGraph*)new TGraphErrors() : (TGraph*)new TGraphAsymmErrors(), ic);
2066 gm->SetLineColor(kBlack);
2067 gm->SetMarkerStyle(7);
2068 gm->SetMarkerColor(kBlack);
2589cf64 2069 gm->SetNameTitle(Form("m_%d_%02d", ig, ic), "mean");
1ee39b3a 2070 }
1ee39b3a 2071 }
2072 }
2073 }
2074
2075/* printf("\n\n\n"); fGraphS->ls();
2076 printf("\n\n\n"); fGraphM->ls();*/
2077
2078
2079 // DEFINE MODELS
2080 // simple gauss
2081 TF1 fg("fGauss", "gaus", -.5, .5);
2082 // Landau for charge resolution
92c40e64 2083 TF1 fch("fClCh", "landau", 0., 1000.);
2084 // Landau for e+- pt resolution
2085 TF1 fpt("fPt", "landau", -0.1, 0.2);
1ee39b3a 2086
2087 //PROCESS EXPERIMENTAL DISTRIBUTIONS
2088 // Charge resolution
2089 //Process3DL(kCharge, 0, &fl);
2090 // Clusters residuals
92d6d80c 2091 Process3D(kCluster, 0, &fg, 1.e4);
81979445 2092 Process3Dlinked(kCluster, 1, &fg);
6558fd69 2093 fNRefFigures = 3;
1ee39b3a 2094 // Tracklet residual/pulls
6558fd69 2095 Process3D(kTrack , 0, &fg, 1.e4);
81979445 2096 Process3Dlinked(kTrack , 1, &fg);
2097 Process3D(kTrack , 2, &fg, 1.e4);
2098 Process3D(kTrack , 3, &fg);
92d6d80c 2099 Process2D(kTrack , 4, &fg, 1.e3);
6558fd69 2100 fNRefFigures = 7;
a310e49b 2101 // TRDin residual/pulls
92d6d80c 2102 Process3D(kTrackIn, 0, &fg, 1.e4);
81979445 2103 Process3Dlinked(kTrackIn, 1, &fg);
2104 Process3D(kTrackIn, 2, &fg, 1.e4);
2105 Process3D(kTrackIn, 3, &fg);
92d6d80c 2106 Process2D(kTrackIn, 4, &fg, 1.e3);
6558fd69 2107 fNRefFigures = 11;
a310e49b 2108 // TRDout residual/pulls
6558fd69 2109 Process3D(kTrackOut, 0, &fg, 1.e3); // scale to fit - see PlotTrackOut
81979445 2110 Process3Dlinked(kTrackOut, 1, &fg);
2111 Process3D(kTrackOut, 2, &fg, 1.e4);
2112 Process3D(kTrackOut, 3, &fg);
a310e49b 2113 Process2D(kTrackOut, 4, &fg, 1.e3);
6558fd69 2114 fNRefFigures = 15;
1ee39b3a 2115
2116 if(!HasMCdata()) return kTRUE;
2117
2118
2119 //PROCESS MC RESIDUAL DISTRIBUTIONS
2120
2121 // CLUSTER Y RESOLUTION/PULLS
92d6d80c 2122 Process3D(kMCcluster, 0, &fg, 1.e4);
81979445 2123 Process3Dlinked(kMCcluster, 1, &fg, 1.);
6558fd69 2124 fNRefFigures = 17;
1ee39b3a 2125
2126 // TRACKLET RESOLUTION/PULLS
92d6d80c 2127 Process3D(kMCtracklet, 0, &fg, 1.e4); // y
81979445 2128 Process3Dlinked(kMCtracklet, 1, &fg, 1.); // y pulls
2129 Process3D(kMCtracklet, 2, &fg, 1.e4); // z
2130 Process3D(kMCtracklet, 3, &fg, 1.); // z pulls
92d6d80c 2131 Process2D(kMCtracklet, 4, &fg, 1.e3); // phi
d25116b6 2132 fNRefFigures = 21;
1ee39b3a 2133
2134 // TRACK RESOLUTION/PULLS
92d6d80c 2135 Process3Darray(kMCtrack, 0, &fg, 1.e4); // y
81979445 2136 Process3DlinkedArray(kMCtrack, 1, &fg); // y PULL
2137 Process3Darray(kMCtrack, 2, &fg, 1.e4); // z
2138 Process3Darray(kMCtrack, 3, &fg); // z PULL
92d6d80c 2139 Process2Darray(kMCtrack, 4, &fg, 1.e3); // phi
2140 Process2Darray(kMCtrack, 5, &fg); // snp PULL
2141 Process2Darray(kMCtrack, 6, &fg, 1.e3); // theta
2142 Process2Darray(kMCtrack, 7, &fg); // tgl PULL
2143 Process3Darray(kMCtrack, 8, &fg, 1.e2); // pt resolution
2144 Process3Darray(kMCtrack, 9, &fg); // 1/pt pulls
2145 Process3Darray(kMCtrack, 10, &fg, 1.e2); // p resolution
d25116b6 2146 fNRefFigures+=16;
92d6d80c 2147
2148 // TRACK TRDin RESOLUTION/PULLS
2149 Process3D(kMCtrackIn, 0, &fg, 1.e4);// y resolution
81979445 2150 Process3Dlinked(kMCtrackIn, 1, &fg); // y pulls
2151 Process3D(kMCtrackIn, 2, &fg, 1.e4);// z resolution
2152 Process3D(kMCtrackIn, 3, &fg); // z pulls
92d6d80c 2153 Process2D(kMCtrackIn, 4, &fg, 1.e3);// phi resolution
2154 Process2D(kMCtrackIn, 5, &fg); // snp pulls
2155 Process2D(kMCtrackIn, 6, &fg, 1.e3);// theta resolution
2156 Process2D(kMCtrackIn, 7, &fg); // tgl pulls
2157 Process3D(kMCtrackIn, 8, &fg, 1.e2);// pt resolution
2158 Process3D(kMCtrackIn, 9, &fg); // 1/pt pulls
2159 Process3D(kMCtrackIn, 10, &fg, 1.e2);// p resolution
d25116b6 2160 fNRefFigures+=8;
92d6d80c 2161
2162 // TRACK TRDout RESOLUTION/PULLS
2163 Process3D(kMCtrackOut, 0, &fg, 1.e4);// y resolution
81979445 2164 Process3Dlinked(kMCtrackOut, 1, &fg); // y pulls
2165 Process3D(kMCtrackOut, 2, &fg, 1.e4);// z resolution
2166 Process3D(kMCtrackOut, 3, &fg); // z pulls
92d6d80c 2167 Process2D(kMCtrackOut, 4, &fg, 1.e3);// phi resolution
2168 Process2D(kMCtrackOut, 5, &fg); // snp pulls
2169 Process2D(kMCtrackOut, 6, &fg, 1.e3);// theta resolution
2170 Process2D(kMCtrackOut, 7, &fg); // tgl pulls
2171 Process3D(kMCtrackOut, 8, &fg, 1.e2);// pt resolution
2172 Process3D(kMCtrackOut, 9, &fg); // 1/pt pulls
2173 Process3D(kMCtrackOut, 10, &fg, 1.e2);// p resolution
d25116b6 2174 fNRefFigures+=8;
1ee39b3a 2175
2176 return kTRUE;
2177}
2178
2179
2180//________________________________________________________
2181void AliTRDresolution::Terminate(Option_t *opt)
2182{
2183 AliTRDrecoTask::Terminate(opt);
2184 if(HasPostProcess()) PostProcess();
2185}
2186
2187//________________________________________________________
2188void AliTRDresolution::AdjustF1(TH1 *h, TF1 *f)
2189{
2190// Helper function to avoid duplication of code
2191// Make first guesses on the fit parameters
2192
2193 // find the intial parameters of the fit !! (thanks George)
2194 Int_t nbinsy = Int_t(.5*h->GetNbinsX());
2195 Double_t sum = 0.;
2196 for(Int_t jbin=nbinsy-4; jbin<=nbinsy+4; jbin++) sum+=h->GetBinContent(jbin); sum/=9.;
2197 f->SetParLimits(0, 0., 3.*sum);
2198 f->SetParameter(0, .9*sum);
2199 Double_t rms = h->GetRMS();
2200 f->SetParLimits(1, -rms, rms);
2201 f->SetParameter(1, h->GetMean());
2202
2203 f->SetParLimits(2, 0., 2.*rms);
2204 f->SetParameter(2, rms);
2205 if(f->GetNpar() <= 4) return;
2206
2207 f->SetParLimits(3, 0., sum);
2208 f->SetParameter(3, .1*sum);
2209
2210 f->SetParLimits(4, -.3, .3);
2211 f->SetParameter(4, 0.);
2212
2213 f->SetParLimits(5, 0., 1.e2);
2214 f->SetParameter(5, 2.e-1);
2215}
2216
afca20ef 2217//________________________________________________________
2589cf64 2218TObjArray* AliTRDresolution::BuildMonitorContainerCluster(const char* name, Bool_t expand)
afca20ef 2219{
3d2a3dff 2220// Build performance histograms for AliTRDcluster.vs TRD track or MC
afca20ef 2221// - y reziduals/pulls
3d2a3dff 2222
2223 TObjArray *arr = new TObjArray(2);
afca20ef 2224 arr->SetName(name); arr->SetOwner();
2225 TH1 *h(NULL); char hname[100], htitle[300];
2226
2227 // tracklet resolution/pull in y direction
3d2a3dff 2228 sprintf(hname, "%s_%s_Y", GetNameId(), name);
2589cf64 2229 sprintf(htitle, "Y res for \"%s\" @ %s;tg(#phi);#Delta y [cm];%s", GetNameId(), name, fgkResYsegmName[fSegmentLevel]);
3d2a3dff 2230 if(!(h = (TH3S*)gROOT->FindObject(hname))){
2589cf64 2231 Int_t nybins=fgkNresYsegm[fSegmentLevel];
2232 if(expand) nybins*=2;
3d2a3dff 2233 h = new TH3S(hname, htitle,
2589cf64 2234 48, -.48, .48, 60, -.15, .15, nybins, -0.5, nybins-0.5);
afca20ef 2235 } else h->Reset();
2236 arr->AddAt(h, 0);
81979445 2237 sprintf(hname, "%s_%s_YZpull", GetNameId(), name);
2589cf64 2238 sprintf(htitle, "YZ pull for \"%s\" @ %s;%s;#Delta y / #sigma_{y};#Delta z / #sigma_{z}", GetNameId(), name, fgkResYsegmName[fSegmentLevel]);
81979445 2239 if(!(h = (TH3S*)gROOT->FindObject(hname))){
2589cf64 2240 h = new TH3S(hname, htitle, fgkNresYsegm[fSegmentLevel], -0.5, fgkNresYsegm[fSegmentLevel]-0.5, 100, -4.5, 4.5, 100, -4.5, 4.5);
afca20ef 2241 } else h->Reset();
2242 arr->AddAt(h, 1);
2243
3d2a3dff 2244 return arr;
2245}
2246
2247//________________________________________________________
2589cf64 2248TObjArray* AliTRDresolution::BuildMonitorContainerTracklet(const char* name, Bool_t expand)
3d2a3dff 2249{
2250// Build performance histograms for AliExternalTrackParam.vs TRD tracklet
2251// - y reziduals/pulls
2252// - z reziduals/pulls
2253// - phi reziduals
2589cf64 2254 TObjArray *arr = BuildMonitorContainerCluster(name, expand);
3d2a3dff 2255 arr->Expand(5);
2256 TH1 *h(NULL); char hname[100], htitle[300];
2257
afca20ef 2258 // tracklet resolution/pull in z direction
3d2a3dff 2259 sprintf(hname, "%s_%s_Z", GetNameId(), name);
81979445 2260 sprintf(htitle, "Z res for \"%s\" @ %s;tg(#theta);#Delta z [cm];row cross", GetNameId(), name);
2261 if(!(h = (TH3S*)gROOT->FindObject(hname))){
2262 h = new TH3S(hname, htitle, 50, -1., 1., 100, -1.5, 1.5, 2, -0.5, 1.5);
afca20ef 2263 } else h->Reset();
2264 arr->AddAt(h, 2);
3d2a3dff 2265 sprintf(hname, "%s_%s_Zpull", GetNameId(), name);
81979445 2266 sprintf(htitle, "Z pull for \"%s\" @ %s;tg(#theta);#Delta z / #sigma_{z};row cross", GetNameId(), name);
2267 if(!(h = (TH3S*)gROOT->FindObject(hname))){
2268 h = new TH3S(hname, htitle, 50, -1., 1., 100, -5.5, 5.5, 2, -0.5, 1.5);
dfd7d48b 2269 h->GetZaxis()->SetBinLabel(1, "no RC");
2270 h->GetZaxis()->SetBinLabel(2, "RC");
afca20ef 2271 } else h->Reset();
2272 arr->AddAt(h, 3);
2273
2274 // tracklet to track phi resolution
3d2a3dff 2275 sprintf(hname, "%s_%s_PHI", GetNameId(), name);
2276 sprintf(htitle, "#Phi res for \"%s\" @ %s;tg(#phi);#Delta #phi [rad];entries", GetNameId(), name);
afca20ef 2277 if(!(h = (TH2I*)gROOT->FindObject(hname))){
2278 h = new TH2I(hname, htitle, 21, -.33, .33, 100, -.5, .5);
2279 } else h->Reset();
2280 arr->AddAt(h, 4);
2281
2282 return arr;
2283}
2284
2285//________________________________________________________
2286TObjArray* AliTRDresolution::BuildMonitorContainerTrack(const char* name)
2287{
2288// Build performance histograms for AliExternalTrackParam.vs MC
2289// - y resolution/pulls
2290// - z resolution/pulls
2291// - phi resolution, snp pulls
2292// - theta resolution, tgl pulls
2293// - pt resolution, 1/pt pulls, p resolution
2294
afca20ef 2295 TObjArray *arr = BuildMonitorContainerTracklet(name);
2296 arr->Expand(11);
3d2a3dff 2297 TH1 *h(NULL); char hname[100], htitle[300];
a310e49b 2298 TAxis *ax(NULL);
3d2a3dff 2299
afca20ef 2300 // snp pulls
3d2a3dff 2301 sprintf(hname, "%s_%s_SNPpull", GetNameId(), name);
2302 sprintf(htitle, "SNP pull for \"%s\" @ %s;tg(#phi);#Delta snp / #sigma_{snp};entries", GetNameId(), name);
afca20ef 2303 if(!(h = (TH2I*)gROOT->FindObject(hname))){
2304 h = new TH2I(hname, htitle, 60, -.3, .3, 100, -4.5, 4.5);
2305 } else h->Reset();
2306 arr->AddAt(h, 5);
2307
2308 // theta resolution
3d2a3dff 2309 sprintf(hname, "%s_%s_THT", GetNameId(), name);
2310 sprintf(htitle, "#Theta res for \"%s\" @ %s;tg(#theta);#Delta #theta [rad];entries", GetNameId(), name);
afca20ef 2311 if(!(h = (TH2I*)gROOT->FindObject(hname))){
2312 h = new TH2I(hname, htitle, 100, -1., 1., 100, -5e-3, 5e-3);
2313 } else h->Reset();
2314 arr->AddAt(h, 6);
2315 // tgl pulls
3d2a3dff 2316 sprintf(hname, "%s_%s_TGLpull", GetNameId(), name);
2317 sprintf(htitle, "TGL pull for \"%s\" @ %s;tg(#theta);#Delta tgl / #sigma_{tgl};entries", GetNameId(), name);
afca20ef 2318 if(!(h = (TH2I*)gROOT->FindObject(hname))){
2319 h = new TH2I(hname, htitle, 100, -1., 1., 100, -4.5, 4.5);
2320 } else h->Reset();
2321 arr->AddAt(h, 7);
2322
2323 const Int_t kNpt(14);
2324 const Int_t kNdpt(150);
2325 const Int_t kNspc = 2*AliPID::kSPECIES+1;
2326 Float_t Pt=0.1, DPt=-.1, Spc=-5.5;
2327 Float_t binsPt[kNpt+1], binsSpc[kNspc+1], binsDPt[kNdpt+1];
2328 for(Int_t i=0;i<kNpt+1; i++,Pt=TMath::Exp(i*.15)-1.) binsPt[i]=Pt;
2329 for(Int_t i=0; i<kNspc+1; i++,Spc+=1.) binsSpc[i]=Spc;
2330 for(Int_t i=0; i<kNdpt+1; i++,DPt+=2.e-3) binsDPt[i]=DPt;
2331
2332 // Pt resolution
3d2a3dff 2333 sprintf(hname, "%s_%s_Pt", GetNameId(), name);
2334 sprintf(htitle, "P_{t} res for \"%s\" @ %s;p_{t} [GeV/c];#Delta p_{t}/p_{t}^{MC};SPECIES", GetNameId(), name);
afca20ef 2335 if(!(h = (TH3S*)gROOT->FindObject(hname))){
2336 h = new TH3S(hname, htitle,
2337 kNpt, binsPt, kNdpt, binsDPt, kNspc, binsSpc);
a310e49b 2338 ax = h->GetZaxis();
2339 for(Int_t ib(1); ib<=ax->GetNbins(); ib++) ax->SetBinLabel(ib, fgParticle[ib-1]);
afca20ef 2340 } else h->Reset();
2341 arr->AddAt(h, 8);
2342 // 1/Pt pulls
3d2a3dff 2343 sprintf(hname, "%s_%s_1Pt", GetNameId(), name);
2344 sprintf(htitle, "1/P_{t} pull for \"%s\" @ %s;1/p_{t}^{MC} [c/GeV];#Delta(1/p_{t})/#sigma(1/p_{t});SPECIES", GetNameId(), name);
afca20ef 2345 if(!(h = (TH3S*)gROOT->FindObject(hname))){
2346 h = new TH3S(hname, htitle,
2347 kNpt, 0., 2., 100, -4., 4., kNspc, -5.5, 5.5);
a310e49b 2348 ax = h->GetZaxis();
2349 for(Int_t ib(1); ib<=ax->GetNbins(); ib++) ax->SetBinLabel(ib, fgParticle[ib-1]);
afca20ef 2350 } else h->Reset();
2351 arr->AddAt(h, 9);
2352 // P resolution
3d2a3dff 2353 sprintf(hname, "%s_%s_P", GetNameId(), name);
2354 sprintf(htitle, "P res for \"%s\" @ %s;p [GeV/c];#Delta p/p^{MC};SPECIES", GetNameId(), name);
afca20ef 2355 if(!(h = (TH3S*)gROOT->FindObject(hname))){
2356 h = new TH3S(hname, htitle,
2357 kNpt, binsPt, kNdpt, binsDPt, kNspc, binsSpc);
a310e49b 2358 ax = h->GetZaxis();
2359 for(Int_t ib(1); ib<=ax->GetNbins(); ib++) ax->SetBinLabel(ib, fgParticle[ib-1]);
afca20ef 2360 } else h->Reset();
2361 arr->AddAt(h, 10);
2362
2363 return arr;
2364}
2365
2366
1ee39b3a 2367//________________________________________________________
2368TObjArray* AliTRDresolution::Histos()
2369{
2370 //
2371 // Define histograms
2372 //
2373
2374 if(fContainer) return fContainer;
2375
92c40e64 2376 fContainer = new TObjArray(kNviews);
1ee39b3a 2377 //fContainer->SetOwner(kTRUE);
31c8fa8a 2378 TH1 *h(NULL);
2379 TObjArray *arr(NULL);
1ee39b3a 2380
31c8fa8a 2381 // binnings for plots containing momentum or pt
2382 const Int_t kNpt(14), kNphi(48), kNdy(60);
2383 Float_t Phi=-.48, Dy=-.3, Pt=0.1;
2384 Float_t binsPhi[kNphi+1], binsDy[kNdy+1], binsPt[kNpt+1];
2385 for(Int_t i=0; i<kNphi+1; i++,Phi+=.02) binsPhi[i]=Phi;
2386 for(Int_t i=0; i<kNdy+1; i++,Dy+=.01) binsDy[i]=Dy;
2387 for(Int_t i=0;i<kNpt+1; i++,Pt=TMath::Exp(i*.15)-1.) binsPt[i]=Pt;
df0514f6 2388
1ee39b3a 2389 // cluster to track residuals/pulls
92d6d80c 2390 fContainer->AddAt(arr = new TObjArray(2), kCharge);
1ee39b3a 2391 arr->SetName("Charge");
2392 if(!(h = (TH3S*)gROOT->FindObject("hCharge"))){
2393 h = new TH3S("hCharge", "Charge Resolution", 20, 1., 2., 24, 0., 3.6, 100, 0., 500.);
2394 h->GetXaxis()->SetTitle("dx/dx_{0}");
2395 h->GetYaxis()->SetTitle("x_{d} [cm]");
2396 h->GetZaxis()->SetTitle("dq/dx [ADC/cm]");
2397 } else h->Reset();
2398 arr->AddAt(h, 0);
2399
2400 // cluster to track residuals/pulls
3d2a3dff 2401 fContainer->AddAt(BuildMonitorContainerCluster("Cl"), kCluster);
afca20ef 2402 // tracklet to TRD track
2589cf64 2403 fContainer->AddAt(BuildMonitorContainerTracklet("Trk", kTRUE), kTrack);
afca20ef 2404 // tracklet to TRDin
2589cf64 2405 fContainer->AddAt(BuildMonitorContainerTracklet("TrkIN", kTRUE), kTrackIn);
a310e49b 2406 // tracklet to TRDout
2407 fContainer->AddAt(BuildMonitorContainerTracklet("TrkOUT"), kTrackOut);
1ee39b3a 2408
2409
2410 // Resolution histos
2411 if(!HasMCdata()) return fContainer;
2412
3d2a3dff 2413 // cluster resolution
2414 fContainer->AddAt(BuildMonitorContainerCluster("MCcl"), kMCcluster);
1ee39b3a 2415
3d2a3dff 2416 // tracklet resolution
2417 fContainer->AddAt(BuildMonitorContainerTracklet("MCtracklet"), kMCtracklet);
1ee39b3a 2418
3d2a3dff 2419 // track resolution
92d6d80c 2420 fContainer->AddAt(arr = new TObjArray(AliTRDgeometry::kNlayer), kMCtrack);
3d2a3dff 2421 arr->SetName("MCtrk");
2422 for(Int_t il(0); il<AliTRDgeometry::kNlayer; il++) arr->AddAt(BuildMonitorContainerTrack(Form("MCtrk_Ly%d", il)), il);
31c8fa8a 2423
afca20ef 2424 // TRDin TRACK RESOLUTION
92d6d80c 2425 fContainer->AddAt(BuildMonitorContainerTrack("MCtrkIN"), kMCtrackIn);
1ee39b3a 2426
afca20ef 2427 // TRDout TRACK RESOLUTION
92d6d80c 2428 fContainer->AddAt(BuildMonitorContainerTrack("MCtrkOUT"), kMCtrackOut);
1ee39b3a 2429
2430 return fContainer;
2431}
2432
2589cf64 2433//________________________________________________________
fe1d1beb 2434Bool_t AliTRDresolution::Load(const Char_t *file, const Char_t *dir)
2589cf64 2435{
2436// Custom load function. Used to guess the segmentation level of the data.
2437
fe1d1beb 2438 if(!AliTRDrecoTask::Load(file, dir)) return kFALSE;
2589cf64 2439
2440 // look for cluster residual plot - always available
2441 TH3S* h3((TH3S*)((TObjArray*)fContainer->At(kClToTrk))->At(0));
2442 Int_t segmentation(h3->GetNbinsZ()/2);
2443 if(segmentation==fgkNresYsegm[0]){ // default segmentation. Nothing to do
2444 return kTRUE;
2445 } else if(segmentation==fgkNresYsegm[1]){ // stack segmentation.
2446 SetSegmentationLevel(1);
2447 } else if(segmentation==fgkNresYsegm[2]){ // detector segmentation.
2448 SetSegmentationLevel(2);
2449 } else {
2450 AliError(Form("Unknown segmentation [%d].", h3->GetNbinsZ()));
2451 return kFALSE;
2452 }
2453
2454 AliDebug(2, Form("Segmentation set to level \"%s\"", fgkResYsegmName[fSegmentLevel]));
2455 return kTRUE;
2456}
2457
2458
1ee39b3a 2459//________________________________________________________
2460Bool_t AliTRDresolution::Process(TH2 * const h2, TF1 *f, Float_t k, TGraphErrors **g)
2461{
2462 //
2463 // Do the processing
2464 //
2465
92c40e64 2466 Char_t pn[10]; sprintf(pn, "p%03d", fIdxPlot);
1ee39b3a 2467 Int_t n = 0;
2468 if((n=g[0]->GetN())) for(;n--;) g[0]->RemovePoint(n);
2469 if((n=g[1]->GetN())) for(;n--;) g[1]->RemovePoint(n);
2589cf64 2470 if(Int_t(h2->GetEntries())){
2471 AliDebug(4, Form("%s: g[%s %s]", pn, g[0]->GetName(), g[0]->GetTitle()));
2472 } else {
2473 AliDebug(2, Form("%s: g[%s %s]: Missing entries.", pn, g[0]->GetName(), g[0]->GetTitle()));
2474 fIdxPlot++;
2475 return kTRUE;
2476 }
1ee39b3a 2477
dfd7d48b 2478 const Int_t kINTEGRAL=1;
2479 for(Int_t ibin = 0; ibin < Int_t(h2->GetNbinsX()/kINTEGRAL); ibin++){
2480 Int_t abin(ibin*kINTEGRAL+1),bbin(abin+kINTEGRAL-1),mbin(abin+Int_t(kINTEGRAL/2));
2481 Double_t x = h2->GetXaxis()->GetBinCenter(mbin);
2482 TH1D *h = h2->ProjectionY(pn, abin, bbin);
068e2c00 2483 if((n=(Int_t)h->GetEntries())<150){
2589cf64 2484 AliDebug(4, Form(" x[%f] range[%d %d] stat[%d] low statistics !", x, abin, bbin, n));
2485 continue;
2486 }
1ee39b3a 2487 h->Fit(f, "QN");
1ee39b3a 2488 Int_t ip = g[0]->GetN();
2589cf64 2489 AliDebug(4, Form(" x_%d[%f] range[%d %d] stat[%d] M[%f] Sgm[%f]", ip, x, abin, bbin, n, f->GetParameter(1), f->GetParameter(2)));
1ee39b3a 2490 g[0]->SetPoint(ip, x, k*f->GetParameter(1));
2491 g[0]->SetPointError(ip, 0., k*f->GetParError(1));
2492 g[1]->SetPoint(ip, x, k*f->GetParameter(2));
2493 g[1]->SetPointError(ip, 0., k*f->GetParError(2));
1ee39b3a 2494/*
2495 g[0]->SetPoint(ip, x, k*h->GetMean());
2496 g[0]->SetPointError(ip, 0., k*h->GetMeanError());
2497 g[1]->SetPoint(ip, x, k*h->GetRMS());
2498 g[1]->SetPointError(ip, 0., k*h->GetRMSError());*/
2499 }
2500 fIdxPlot++;
2501 return kTRUE;
2502}
2503
2504//________________________________________________________
92c40e64 2505Bool_t AliTRDresolution::Process2D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k, Int_t gidx)
1ee39b3a 2506{
2507 //
2508 // Do the processing
2509 //
2510
2511 if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
2512
2513 // retrive containers
a310e49b 2514 TH2I *h2(NULL);
2515 if(idx<0){
2516 if(!(h2= (TH2I*)(fContainer->At(plot)))) return kFALSE;
2517 } else{
2518 TObjArray *a0(NULL);
2519 if(!(a0=(TObjArray*)(fContainer->At(plot)))) return kFALSE;
2520 if(!(h2=(TH2I*)a0->At(idx))) return kFALSE;
2521 }
2589cf64 2522 if(Int_t(h2->GetEntries())){
2523 AliDebug(2, Form("p[%d] idx[%d] : h[%s] %s", plot, idx, h2->GetName(), h2->GetTitle()));
2524 } else {
2525 AliDebug(2, Form("p[%d] idx[%d] : Missing entries.", plot, idx));
2526 return kFALSE;
2527 }
92c40e64 2528
1ee39b3a 2529 TGraphErrors *g[2];
92c40e64 2530 if(gidx<0) gidx=idx;
2531 if(!(g[0] = gidx<0 ? (TGraphErrors*)fGraphM->At(plot) : (TGraphErrors*)((TObjArray*)(fGraphM->At(plot)))->At(gidx))) return kFALSE;
1ee39b3a 2532
92c40e64 2533 if(!(g[1] = gidx<0 ? (TGraphErrors*)fGraphS->At(plot) : (TGraphErrors*)((TObjArray*)(fGraphS->At(plot)))->At(gidx))) return kFALSE;
1ee39b3a 2534
2535 return Process(h2, f, k, g);
2536}
2537
2538//________________________________________________________
2539Bool_t AliTRDresolution::Process3D(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k)
2540{
2541 //
2542 // Do the processing
2543 //
2544
2545 if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
2546
2547 // retrive containers
a310e49b 2548 TH3S *h3(NULL);
2549 if(idx<0){
2550 if(!(h3= (TH3S*)(fContainer->At(plot)))) return kFALSE;
2551 } else{
2552 TObjArray *a0(NULL);
2553 if(!(a0=(TObjArray*)(fContainer->At(plot)))) return kFALSE;
2554 if(!(h3=(TH3S*)a0->At(idx))) return kFALSE;
2555 }
2589cf64 2556 if(Int_t(h3->GetEntries())){
2557 AliDebug(2, Form("p[%d] idx[%d] h[%s] %s", plot, idx, h3->GetName(), h3->GetTitle()));
2558 } else {
2559 AliDebug(2, Form("p[%d] idx[%d] : Missing entries.", plot, idx));
2560 return kFALSE;
2561 }
1ee39b3a 2562
2563 TObjArray *gm, *gs;
2564 if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE;
2565 if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE;
2566 TGraphErrors *g[2];
2567
2568 TAxis *az = h3->GetZaxis();
dfd7d48b 2569 for(Int_t iz(0); iz<gm->GetEntriesFast(); iz++){
2570 if(!(g[0] = (TGraphErrors*)gm->At(iz))) return kFALSE;
2571 if(!(g[1] = (TGraphErrors*)gs->At(iz))) return kFALSE;
2572 az->SetRange(iz+1, iz+1);
1ee39b3a 2573 if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE;
2574 }
2575
2576 return kTRUE;
2577}
2578
92c40e64 2579
81979445 2580//________________________________________________________
2581Bool_t AliTRDresolution::Process3Dlinked(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k)
2582{
2583 //
2584 // Do the processing
2585 //
2586
2587 if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
2588
2589 // retrive containers
2590 TH3S *h3(NULL);
2591 if(idx<0){
2592 if(!(h3= (TH3S*)(fContainer->At(plot)))) return kFALSE;
2593 } else{
2594 TObjArray *a0(NULL);
2595 if(!(a0=(TObjArray*)(fContainer->At(plot)))) return kFALSE;
2596 if(!(h3=(TH3S*)a0->At(idx))) return kFALSE;
2597 }
2589cf64 2598 if(Int_t(h3->GetEntries())){
2599 AliDebug(2, Form("p[%d] idx[%d] h[%s] %s", plot, idx, h3->GetName(), h3->GetTitle()));
2600 } else {
2601 AliDebug(2, Form("p[%d] idx[%d] : Missing entries.", plot, idx));
2602 return kFALSE;
2603 }
81979445 2604
2605 TObjArray *gm, *gs;
2606 if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE;
2607 if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE;
2608 TGraphErrors *g[2];
2609
2610 if(!(g[0] = (TGraphErrors*)gm->At(0))) return kFALSE;
2611 if(!(g[1] = (TGraphErrors*)gs->At(0))) return kFALSE;
2612 if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE;
2613
2614 if(!(g[0] = (TGraphErrors*)gm->At(1))) return kFALSE;
2615 if(!(g[1] = (TGraphErrors*)gs->At(1))) return kFALSE;
2616 if(!Process((TH2*)h3->Project3D("zx"), f, k, g)) return kFALSE;
2617
2618 return kTRUE;
2619}
2620
2621
1ee39b3a 2622//________________________________________________________
2623Bool_t AliTRDresolution::Process3DL(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k)
2624{
2625 //
2626 // Do the processing
2627 //
2628
2629 if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
2630
2631 // retrive containers
2632 TH3S *h3 = (TH3S*)((TObjArray*)fContainer->At(plot))->At(idx);
2633 if(!h3) return kFALSE;
92d6d80c 2634 AliDebug(2, Form("p[%d] idx[%d] h[%s] %s", plot, idx, h3->GetName(), h3->GetTitle()));
1ee39b3a 2635
2636 TGraphAsymmErrors *gm;
2637 TGraphErrors *gs;
2638 if(!(gm = (TGraphAsymmErrors*)((TObjArray*)fGraphM->At(plot))->At(0))) return kFALSE;
2639 if(!(gs = (TGraphErrors*)((TObjArray*)fGraphS->At(plot)))) return kFALSE;
2640
2641 Float_t x, r, mpv, xM, xm;
2642 TAxis *ay = h3->GetYaxis();
2643 for(Int_t iy=1; iy<=h3->GetNbinsY(); iy++){
2644 ay->SetRange(iy, iy);
2645 x = ay->GetBinCenter(iy);
2646 TH2F *h2=(TH2F*)h3->Project3D("zx");
2647 TAxis *ax = h2->GetXaxis();
2648 for(Int_t ix=1; ix<=h2->GetNbinsX(); ix++){
2649 r = ax->GetBinCenter(ix);
2650 TH1D *h1 = h2->ProjectionY("py", ix, ix);
2651 if(h1->Integral()<50) continue;
2652 h1->Fit(f, "QN");
2653
2654 GetLandauMpvFwhm(f, mpv, xm, xM);
2655 Int_t ip = gm->GetN();
2656 gm->SetPoint(ip, x, k*mpv);
2657 gm->SetPointError(ip, 0., 0., k*xm, k*xM);
2658 gs->SetPoint(ip, r, k*(xM-xm)/mpv);
2659 gs->SetPointError(ip, 0., 0.);
2660 }
2661 }
2662
2663 return kTRUE;
2664}
2665
2666//________________________________________________________
92d6d80c 2667Bool_t AliTRDresolution::Process2Darray(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k)
2668{
2669 //
2670 // Do the processing
2671 //
2672
2673 if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
2674
2675 // retrive containers
2676 TObjArray *arr = (TObjArray*)(fContainer->At(plot));
2677 if(!arr) return kFALSE;
2678 AliDebug(2, Form("p[%d] idx[%d] arr[%s]", plot, idx, arr->GetName()));
2679
2680 TObjArray *gm, *gs;
2681 if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE;
2682 if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE;
2683
2684 TGraphErrors *g[2]; TH2I *h2(NULL); TObjArray *a0(NULL);
2685 for(Int_t ia(0); ia<arr->GetEntriesFast(); ia++){
2686 if(!(a0 = (TObjArray*)arr->At(ia))) continue;
2687
2688 if(!(h2 = (TH2I*)a0->At(idx))) return kFALSE;
2589cf64 2689 if(Int_t(h2->GetEntries())){
2690 AliDebug(4, Form(" idx[%d] h[%s] %s", ia, h2->GetName(), h2->GetTitle()));
2691 } else {
2692 AliDebug(2, Form(" idx[%d] : Missing entries.", ia));
2693 continue;
2694 }
92d6d80c 2695
2696 if(!(g[0] = (TGraphErrors*)gm->At(ia))) return kFALSE;
2697 if(!(g[1] = (TGraphErrors*)gs->At(ia))) return kFALSE;
2698 if(!Process(h2, f, k, g)) return kFALSE;
2699 }
2700
2701 return kTRUE;
2702}
2703
2704//________________________________________________________
2705Bool_t AliTRDresolution::Process3Darray(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k)
1ee39b3a 2706{
2707 //
2708 // Do the processing
2709 //
2710
2711 if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
92c40e64 2712 //printf("Process4D : processing plot[%d] idx[%d]\n", plot, idx);
1ee39b3a 2713
2714 // retrive containers
92d6d80c 2715 TObjArray *arr = (TObjArray*)(fContainer->At(plot));
1ee39b3a 2716 if(!arr) return kFALSE;
92d6d80c 2717 AliDebug(2, Form("p[%d] idx[%d] arr[%s]", plot, idx, arr->GetName()));
1ee39b3a 2718
92c40e64 2719 TObjArray *gm, *gs;
2720 if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE;
2721 if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE;
1ee39b3a 2722
92d6d80c 2723 TGraphErrors *g[2]; TH3S *h3(NULL); TObjArray *a0(NULL);
2724 Int_t in(0);
2725 for(Int_t ia(0); ia<arr->GetEntriesFast(); ia++){
2726 if(!(a0 = (TObjArray*)arr->At(ia))) continue;
1ee39b3a 2727
92d6d80c 2728 if(!(h3 = (TH3S*)a0->At(idx))) return kFALSE;
2589cf64 2729 if(Int_t(h3->GetEntries())){
2730 AliDebug(4, Form(" idx[%d] h[%s] %s", ia, h3->GetName(), h3->GetTitle()));
2731 } else {
2732 AliDebug(2, Form(" idx[%d] : Missing entries.", ia));
2733 continue;
2734 }
1ee39b3a 2735 TAxis *az = h3->GetZaxis();
92c40e64 2736 for(Int_t iz=1; iz<=az->GetNbins(); iz++, in++){
dfd7d48b 2737 if(in >= gm->GetEntriesFast()) break;
92c40e64 2738 if(!(g[0] = (TGraphErrors*)gm->At(in))) return kFALSE;
2739 if(!(g[1] = (TGraphErrors*)gs->At(in))) return kFALSE;
1ee39b3a 2740 az->SetRange(iz, iz);
2741 if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE;
2742 }
2743 }
92d6d80c 2744 AliDebug(2, Form("Projections [%d] from [%d]", in, gs->GetEntriesFast()));
1ee39b3a 2745
2746 return kTRUE;
2747}
2748
81979445 2749//________________________________________________________
2750Bool_t AliTRDresolution::Process3DlinkedArray(ETRDresolutionPlot plot, Int_t idx, TF1 *f, Float_t k)
2751{
2752 //
2753 // Do the processing
2754 //
2755
2756 if(!fContainer || !fGraphS || !fGraphM) return kFALSE;
2757 //printf("Process4D : processing plot[%d] idx[%d]\n", plot, idx);
2758
2759 // retrive containers
2760 TObjArray *arr = (TObjArray*)(fContainer->At(plot));
2761 if(!arr) return kFALSE;
2762 AliDebug(2, Form("p[%d] idx[%d] arr[%s]", plot, idx, arr->GetName()));
2763
2764 TObjArray *gm, *gs;
2765 if(!(gm = (TObjArray*)((TObjArray*)(fGraphM->At(plot)))->At(idx))) return kFALSE;
2766 if(!(gs = (TObjArray*)((TObjArray*)(fGraphS->At(plot)))->At(idx))) return kFALSE;
2767
2768 TGraphErrors *g[2]; TH3S *h3(NULL); TObjArray *a0(NULL);
2769 Int_t in(0);
2770 for(Int_t ia(0); ia<arr->GetEntriesFast(); ia++){
2771 if(!(a0 = (TObjArray*)arr->At(ia))) continue;
2772 if(!(h3 = (TH3S*)a0->At(idx))) return kFALSE;
2589cf64 2773 if(Int_t(h3->GetEntries())){
2774 AliDebug(4, Form(" idx[%d] h[%s] %s", ia, h3->GetName(), h3->GetTitle()));
2775 } else {
2776 AliDebug(2, Form(" idx[%d] : Missing entries.", ia));
2777 continue;
2778 }
81979445 2779 if(!(g[0] = (TGraphErrors*)gm->At(in))) return kFALSE;
2780 if(!(g[1] = (TGraphErrors*)gs->At(in))) return kFALSE;
2781 if(!Process((TH2*)h3->Project3D("yx"), f, k, g)) return kFALSE;
2782 in++;
2783
2784 if(!(g[0] = (TGraphErrors*)gm->At(in))) return kFALSE;
2785 if(!(g[1] = (TGraphErrors*)gs->At(in))) return kFALSE;
2786 if(!Process((TH2*)h3->Project3D("zx"), f, k, g)) return kFALSE;
2787 in++;
2788 }
2789 AliDebug(2, Form("Projections [%d] from [%d]", in, gs->GetEntriesFast()));
2790
2791 return kTRUE;
2792}
2793
1ee39b3a 2794//________________________________________________________
a310e49b 2795Bool_t AliTRDresolution::GetGraph(Float_t *bb, ETRDresolutionPlot ip, Int_t idx, Bool_t kLEG, const Char_t *explain)
1ee39b3a 2796{
2797 //
2798 // Get the graphs
2799 //
2800
2801 if(!fGraphS || !fGraphM) return kFALSE;
a310e49b 2802 // axis titles look up
1ee39b3a 2803 Int_t nref = 0;
92c40e64 2804 for(Int_t jp=0; jp<(Int_t)ip; jp++) nref+=fgNproj[jp];
1ee39b3a 2805 UChar_t jdx = idx<0?0:idx;
92c40e64 2806 for(Int_t jc=0; jc<TMath::Min(jdx,fgNproj[ip]-1); jc++) nref++;
2589cf64 2807 Char_t **at = fAxTitle[nref];
1ee39b3a 2808
a310e49b 2809 // build legends if requiered
2810 TLegend *leg(NULL);
2811 if(kLEG){
2812 leg=new TLegend(.65, .6, .95, .9);
2813 leg->SetBorderSize(0);
2814 leg->SetFillStyle(0);
2815 }
2816 // build frame
2817 TH1S *h1(NULL);
2818 h1 = new TH1S(Form("h1TF_%02d", fIdxFrame++), Form("%s %s;%s;%s", at[0], explain?explain:"", at[1], at[2]), 2, bb[0], bb[2]);
2819 h1->SetMinimum(bb[1]);h1->SetMaximum(bb[3]);
2820 h1->SetLineColor(kBlack); h1->SetLineWidth(1);h1->SetLineStyle(2);
2821 // axis range
2822 TAxis *ax = h1->GetXaxis();
2823 ax->CenterTitle();ax->SetMoreLogLabels();ax->SetTitleOffset(1.2);
2824 ax = h1->GetYaxis();
2825 ax->SetRangeUser(bb[1], bb[3]);
2826 ax->CenterTitle(); ax->SetTitleOffset(1.4);
2827 h1->Draw();
2828 // bounding box
2829 TBox *b = new TBox(-.15, bb[1], .15, bb[3]);
2830 b->SetFillStyle(3002);b->SetLineColor(0);
2831 b->SetFillColor(ip<=Int_t(kMCcluster)?kGreen:kBlue);
2832 b->Draw();
2833
2834 TGraphErrors *gm = idx<0 ? (TGraphErrors*)fGraphM->At(ip) : (TGraphErrors*)((TObjArray*)(fGraphM->At(ip)))->At(idx);
2835 if(!gm) return kFALSE;
2836 TGraphErrors *gs = idx<0 ? (TGraphErrors*)fGraphS->At(ip) : (TGraphErrors*)((TObjArray*)(fGraphS->At(ip)))->At(idx);
2837 if(!gs) return kFALSE;
2838
2839 Int_t n(0), nPlots(0);
1ee39b3a 2840 if((n=gm->GetN())) {
a310e49b 2841 nPlots++;
2842 gm->Draw("pl"); if(leg) leg->AddEntry(gm, gm->GetTitle(), "pl");
1ee39b3a 2843 PutTrendValue(Form("%s_%s", fgPerformanceName[ip], at[0]), gm->GetMean(2));
2844 PutTrendValue(Form("%s_%sRMS", fgPerformanceName[ip], at[0]), gm->GetRMS(2));
2845 }
2846
2847 if((n=gs->GetN())){
a310e49b 2848 nPlots++;
2849 gs->Draw("pl"); if(leg) leg->AddEntry(gs, gs->GetTitle(), "pl");
1ee39b3a 2850 gs->Sort(&TGraph::CompareY);
2851 PutTrendValue(Form("%s_%sSigMin", fgPerformanceName[ip], at[0]), gs->GetY()[0]);
2852 PutTrendValue(Form("%s_%sSigMax", fgPerformanceName[ip], at[0]), gs->GetY()[n-1]);
2853 gs->Sort(&TGraph::CompareX);
2854 }
a310e49b 2855 if(!nPlots) return kFALSE;
2856 if(leg) leg->Draw();
1ee39b3a 2857 return kTRUE;
2858}
2859
1ee39b3a 2860//________________________________________________________
92d6d80c 2861Bool_t AliTRDresolution::GetGraphArray(Float_t *bb, ETRDresolutionPlot ip, Int_t idx, Bool_t kLEG, Int_t n, Int_t *sel, const Char_t *explain)
1ee39b3a 2862{
2863 //
2864 // Get the graphs
2865 //
2866
2867 if(!fGraphS || !fGraphM) return kFALSE;
2868
2869 // axis titles look up
92c40e64 2870 Int_t nref(0);
92d6d80c 2871 for(Int_t jp(0); jp<ip; jp++) nref+=fgNproj[jp];
92c40e64 2872 nref+=idx;
2589cf64 2873 Char_t **at = fAxTitle[nref];
92c40e64 2874
92d6d80c 2875 // build legends if requiered
31c8fa8a 2876 TLegend *legM(NULL), *legS(NULL);
92c40e64 2877 if(kLEG){
1295b37f 2878 legM=new TLegend(.35, .6, .65, .9);
31c8fa8a 2879 legM->SetHeader("Mean");
92d6d80c 2880 legM->SetBorderSize(0);
2881 legM->SetFillStyle(0);
1295b37f 2882 legS=new TLegend(.65, .6, .95, .9);
31c8fa8a 2883 legS->SetHeader("Sigma");
92d6d80c 2884 legS->SetBorderSize(0);
2885 legS->SetFillStyle(0);
1ee39b3a 2886 }
92d6d80c 2887 // build frame
92c40e64 2888 TH1S *h1(NULL);
92d6d80c 2889 h1 = new TH1S(Form("h1TF_%02d", fIdxFrame++), Form("%s %s;%s;%s", at[0], explain?explain:"", at[1], at[2]), 2, bb[0], bb[2]);
92c40e64 2890 h1->SetMinimum(bb[1]);h1->SetMaximum(bb[3]);
2891 h1->SetLineColor(kBlack); h1->SetLineWidth(1);h1->SetLineStyle(2);
1ee39b3a 2892 // axis range
92c40e64 2893 TAxis *ax = h1->GetXaxis();
1295b37f 2894 ax->CenterTitle();ax->SetMoreLogLabels();ax->SetTitleOffset(1.2);
92c40e64 2895 ax = h1->GetYaxis();
1ee39b3a 2896 ax->SetRangeUser(bb[1], bb[3]);
1295b37f 2897 ax->CenterTitle(); ax->SetTitleOffset(1.4);
92c40e64 2898 h1->Draw();
92d6d80c 2899
2900 TGraphErrors *gm(NULL), *gs(NULL);
2901 TObjArray *a0(NULL), *a1(NULL);
2902 a0 = (TObjArray*)((TObjArray*)fGraphM->At(ip))->At(idx);
2903 a1 = (TObjArray*)((TObjArray*)fGraphS->At(ip))->At(idx);
2904 if(!n) n=a0->GetEntriesFast();
2905 AliDebug(4, Form("Graph : Ref[%d] Title[%s] Limits{x[%f %f] y[%f %f]} Comp[%d] Selection[%c]", nref, at[0], bb[0], bb[2], bb[1], bb[3], n, sel ? 'y' : 'n'));
2906 Int_t nn(0), nPlots(0);
2907 for(Int_t is(0), is0(0); is<n; is++){
2908 is0 = sel ? sel[is] : is;
2909 if(!(gs = (TGraphErrors*)a1->At(is0))) return kFALSE;
2910 if(!(gm = (TGraphErrors*)a0->At(is0))) return kFALSE;
1ee39b3a 2911
92c40e64 2912 if((nn=gs->GetN())){
92d6d80c 2913 nPlots++;
a310e49b 2914 gs->Draw("pc");
2915 if(legS){
2916 //printf("LegEntry %s [%s]%s\n", at[0], gs->GetName(), gs->GetTitle());
2917 legS->AddEntry(gs, gs->GetTitle(), "pl");
2918 }
92c40e64 2919 gs->Sort(&TGraph::CompareY);
92d6d80c 2920 PutTrendValue(Form("%s_%sSigMin", fgPerformanceName[kMCtrack], at[0]), gs->GetY()[0]);
2921 PutTrendValue(Form("%s_%sSigMax", fgPerformanceName[kMCtrack], at[0]), gs->GetY()[nn-1]);
92c40e64 2922 gs->Sort(&TGraph::CompareX);
2923 }
2924 if(gm->GetN()){
92d6d80c 2925 nPlots++;
a310e49b 2926 gm->Draw("pc");
2927 if(legM){
2928 //printf("LegEntry %s [%s]%s\n", at[0], gm->GetName(), gm->GetTitle());
2929 legM->AddEntry(gm, gm->GetTitle(), "pl");
2930 }
92d6d80c 2931 PutTrendValue(Form("%s_%s", fgPerformanceName[kMCtrack], at[0]), gm->GetMean(2));
2932 PutTrendValue(Form("%s_%sRMS", fgPerformanceName[kMCtrack], at[0]), gm->GetRMS(2));
92c40e64 2933 }
2934 }
92d6d80c 2935 if(!nPlots) return kFALSE;
92c40e64 2936 if(kLEG){
92d6d80c 2937 legM->Draw();
2938 legS->Draw();
1ee39b3a 2939 }
1ee39b3a 2940 return kTRUE;
2941}
2942
2943//________________________________________________________
2944void AliTRDresolution::GetLandauMpvFwhm(TF1 * const f, Float_t &mpv, Float_t &xm, Float_t &xM)
2945{
2946 //
2947 // Get the most probable value and the full width half mean
2948 // of a Landau distribution
2949 //
2950
2951 const Float_t dx = 1.;
2952 mpv = f->GetParameter(1);
2953 Float_t fx, max = f->Eval(mpv);
2954
2955 xm = mpv - dx;
2956 while((fx = f->Eval(xm))>.5*max){
2957 if(fx>max){
2958 max = fx;
2959 mpv = xm;
2960 }
2961 xm -= dx;
2962 }
2963
2964 xM += 2*(mpv - xm);
2965 while((fx = f->Eval(xM))>.5*max) xM += dx;
2966}
2967
2968
2969//________________________________________________________
2970void AliTRDresolution::SetRecoParam(AliTRDrecoParam *r)
2971{
2972
2973 fReconstructor->SetRecoParam(r);
2974}
2589cf64 2975
2976
2977//________________________________________________________
2978void AliTRDresolution::SetSegmentationLevel(Int_t l)
2979{
2980// Setting the segmentation level to "l"
2981 fSegmentLevel = l;
2982
2983 UShort_t const lNcomp[kNprojs] = {
2984 1, 1, //2,
2985 fgkNresYsegm[fSegmentLevel], 2, //2,
2986 2*fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5,
2987 2*fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5,
2988 2*fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5,
2989 // MC
2990 fgkNresYsegm[fSegmentLevel], 2, //2,
2991 fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, //5,
2992 fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, 1, 1, 1, 11, 11, 11, //11
2993 fgkNresYsegm[fSegmentLevel], 2, 2, 2, 1, 1, 1, 1, 11, 11, 11, //11
2994 6*fgkNresYsegm[fSegmentLevel], 6*2, 6*2, 6*2, 6, 6, 6, 6, 6*11, 6*11, 6*11 //11
2995 };
2996 memcpy(fNcomp, lNcomp, kNprojs*sizeof(UShort_t));
2997
2998 Char_t const *lAxTitle[kNprojs][4] = {
2999 // Charge
3000 {"Impv", "x [cm]", "I_{mpv}", "x/x_{0}"}
3001 ,{"dI/Impv", "x/x_{0}", "#delta I/I_{mpv}", "x[cm]"}
3002 // Clusters to Kalman
3003 ,{"Cluster2Track residuals", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
3004 ,{"Cluster2Track YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
3005 // TRD tracklet to Kalman fit
3006 ,{"Tracklet2Track Y residuals", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
3007 ,{"Tracklet2Track YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
3008 ,{"Tracklet2Track Z residuals", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
3009 ,{"Tracklet2Track Z pulls", "tg(#theta)", "z", "#sigma_{z}"}
3010 ,{"Tracklet2Track Phi residuals", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
3011 // TRDin 2 first TRD tracklet
3012 ,{"Tracklet2Track Y residuals @ TRDin", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
3013 ,{"Tracklet2Track YZ pulls @ TRDin", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
3014 ,{"Tracklet2Track Z residuals @ TRDin", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
3015 ,{"Tracklet2Track Z pulls @ TRDin", "tg(#theta)", "z", "#sigma_{z}"}
3016 ,{"Tracklet2Track Phi residuals @ TRDin", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
3017 // TRDout 2 first TRD tracklet
3018 ,{"Tracklet2Track Y residuals @ TRDout", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
3019 ,{"Tracklet2Track YZ pulls @ TRDout", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
3020 ,{"Tracklet2Track Z residuals @ TRDout", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
3021 ,{"Tracklet2Track Z pulls @ TRDout", "tg(#theta)", "z", "#sigma_{z}"}
3022 ,{"Tracklet2Track Phi residuals @ TRDout", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
3023 // MC cluster
3024 ,{"MC Cluster Y resolution", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
3025 ,{"MC Cluster YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
3026 // MC tracklet
3027 ,{"MC Tracklet Y resolution", "tg(#phi)", "y [#mum]", "#sigma_{y}[#mum]"}
3028 ,{"MC Tracklet YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
3029 ,{"MC Tracklet Z resolution", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
3030 ,{"MC Tracklet Z pulls", "tg(#theta)", "z", "#sigma_{z}"}
3031 ,{"MC Tracklet Phi resolution", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
3032 // MC track TRDin
3033 ,{"MC Y resolution @ TRDin", "tg(#phi)", "y [#mum]", "#sigma_{y}[#mum]"}
3034 ,{"MC YZ pulls @ TRDin", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
3035 ,{"MC Z resolution @ TRDin", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
3036 ,{"MC Z pulls @ TRDin", "tg(#theta)", "z", "#sigma_{z}"}
3037 ,{"MC #Phi resolution @ TRDin", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
3038 ,{"MC SNP pulls @ TRDin", "tg(#phi)", "SNP", "#sigma_{snp}"}
3039 ,{"MC #Theta resolution @ TRDin", "tg(#theta)", "#theta [mrad]", "#sigma_{#theta} [mrad]"}
3040 ,{"MC TGL pulls @ TRDin", "tg(#theta)", "TGL", "#sigma_{tgl}"}
3041 ,{"MC P_{t} resolution @ TRDin", "p_{t}^{MC} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "MC: #sigma^{TPC}(#Deltap_{t}/p_{t}^{MC}) [%]"}
3042 ,{"MC 1/P_{t} pulls @ TRDin", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC}-1/p_{t}^{MC}", "MC PULL: #sigma_{1/p_{t}}^{TPC}"}
3043 ,{"MC P resolution @ TRDin", "p^{MC} [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "MC: #sigma^{TPC}(#Deltap/p^{MC}) [%]"}
3044 // MC track TRDout
3045 ,{"MC Y resolution @ TRDout", "tg(#phi)", "y [#mum]", "#sigma_{y}[#mum]"}
3046 ,{"MC YZ pulls @ TRDout", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
3047 ,{"MC Z resolution @ TRDout", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
3048 ,{"MC Z pulls @ TRDout", "tg(#theta)", "z", "#sigma_{z}"}
3049 ,{"MC #Phi resolution @ TRDout", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
3050 ,{"MC SNP pulls @ TRDout", "tg(#phi)", "SNP", "#sigma_{snp}"}
3051 ,{"MC #Theta resolution @ TRDout", "tg(#theta)", "#theta [mrad]", "#sigma_{#theta} [mrad]"}
3052 ,{"MC TGL pulls @ TRDout", "tg(#theta)", "TGL", "#sigma_{tgl}"}
3053 ,{"MC P_{t} resolution @ TRDout", "p_{t}^{MC} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "MC: #sigma^{TPC}(#Deltap_{t}/p_{t}^{MC}) [%]"}
3054 ,{"MC 1/P_{t} pulls @ TRDout", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC}-1/p_{t}^{MC}", "MC PULL: #sigma_{1/p_{t}}^{TPC}"}
3055 ,{"MC P resolution @ TRDout", "p^{MC} [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "MC: #sigma^{TPC}(#Deltap/p^{MC}) [%]"}
3056 // MC track in TRD
3057 ,{"MC Track Y resolution", "tg(#phi)", "y [#mum]", "#sigma_{y} [#mum]"}
3058 ,{"MC Track YZ pulls", fgkResYsegmName[fSegmentLevel], "y / z", "#sigma_{y}"}
3059 ,{"MC Track Z resolution", "tg(#theta)", "z [#mum]", "#sigma_{z} [#mum]"}
3060 ,{"MC Track Z pulls", "tg(#theta)", "z", "#sigma_{z}"}
3061 ,{"MC Track #Phi resolution", "tg(#phi)", "#phi [mrad]", "#sigma_{#phi} [mrad]"}
3062 ,{"MC Track SNP pulls", "tg(#phi)", "SNP", "#sigma_{snp}"}
3063 ,{"MC Track #Theta resolution", "tg(#theta)", "#theta [mrad]", "#sigma_{#theta} [mrad]"}
3064 ,{"MC Track TGL pulls", "tg(#theta)", "TGL", "#sigma_{tgl}"}
3065 ,{"MC P_{t} resolution", "p_{t} [GeV/c]", "(p_{t}^{REC}-p_{t}^{MC})/p_{t}^{MC} [%]", "#sigma(#Deltap_{t}/p_{t}^{MC}) [%]"}
3066 ,{"MC 1/P_{t} pulls", "1/p_{t}^{MC} [c/GeV]", "1/p_{t}^{REC} - 1/p_{t}^{MC}", "#sigma_{1/p_{t}}"}
3067 ,{"MC P resolution", "p [GeV/c]", "(p^{REC}-p^{MC})/p^{MC} [%]", "#sigma(#Deltap/p^{MC}) [%]"}
3068 };
3069 memcpy(fAxTitle, lAxTitle, 4*kNprojs*sizeof(Char_t*));
3070}