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