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