]>
Commit | Line | Data |
---|---|---|
e8541476 | 1 | //$Id$ |
2 | ||
b661165c | 3 | // Author: Anders Vestbo <mailto:vestbo@fi.uib.no> |
4 | //*-- Copyright © ASV | |
108615fc | 5 | |
118c26c3 | 6 | #include "AliL3StandardIncludes.h" |
108615fc | 7 | #include <TFile.h> |
8 | #include <TH1.h> | |
9 | #include <TParticle.h> | |
10 | #include <TTree.h> | |
11 | #include <TClonesArray.h> | |
12 | ||
4499ed26 | 13 | #include <AliRun.h> |
14 | #include <AliSimDigits.h> | |
15 | #include <AliTPC.h> | |
16 | #include <AliTPCcluster.h> | |
17 | #include <AliTPCClustersArray.h> | |
18 | #include <AliTPCClustersRow.h> | |
19 | #include <AliTPCParam.h> | |
20 | #include <AliComplexCluster.h> | |
108615fc | 21 | |
118c26c3 | 22 | #if GCCVERSION == 3 |
23 | #include <fstream> | |
24 | #include <iosfwd> | |
25 | #else | |
26 | #include <fstream.h> | |
27 | #endif | |
28 | ||
29 | #include "AliL3Logging.h" | |
108615fc | 30 | #include "AliL3Transform.h" |
31 | #include "AliL3SpacePointData.h" | |
32 | #include "AliL3Track.h" | |
33 | #include "AliL3FileHandler.h" | |
34 | #include "AliL3TrackArray.h" | |
35 | #include "AliL3Evaluate.h" | |
118c26c3 | 36 | |
37 | #if GCCVERSION == 3 | |
38 | using namespace std; | |
39 | #endif | |
108615fc | 40 | |
b661165c | 41 | //_____________________________________________________________ |
42 | // AliL3Evaluate | |
43 | // | |
44 | // Evaluation class for tracking. Plots efficiencies etc.. | |
84122e15 | 45 | // |
46 | ||
108615fc | 47 | |
48 | ClassImp(AliL3Evaluate) | |
49 | ||
50 | AliL3Evaluate::AliL3Evaluate() | |
51 | { | |
e1a32fa5 | 52 | fDigitsFile = NULL; |
108615fc | 53 | fTracks = NULL; |
54 | fMCclusterfile = NULL; | |
c0217812 | 55 | fNFastPoints = 0; |
56 | fMcIndex = 0; | |
57 | fMcId = 0; | |
1f79afc0 | 58 | fMinSlice=0; |
59 | fMaxSlice=0; | |
108615fc | 60 | } |
61 | ||
40896357 | 62 | AliL3Evaluate::AliL3Evaluate(Char_t *mcfile,Int_t *slice) |
108615fc | 63 | { |
64 | //Normal constructor. Input are the rootfile containing the | |
65 | //original MC information of the simulated event. | |
66 | ||
e1a32fa5 | 67 | fEventFile = new TFile(mcfile,"READ"); |
40896357 | 68 | |
5e0f9911 | 69 | fParam = (AliTPCParam*)fEventFile->Get(AliL3Transform::GetParamName()); |
108615fc | 70 | |
71 | fMinSlice = slice[0]; | |
72 | fMaxSlice = slice[1]; | |
e170e0f9 | 73 | fMinGoodPt = 0.1; |
23908b9b | 74 | fNoOverlap = kFALSE; |
108615fc | 75 | } |
76 | ||
1f79afc0 | 77 | AliL3Evaluate::AliL3Evaluate(Int_t *slice) |
78 | { | |
78127e35 | 79 | //ctor to use if you do not need any rootfile. |
80 | ||
1f79afc0 | 81 | fMinSlice = slice[0]; |
82 | fMaxSlice = slice[1]; | |
1f79afc0 | 83 | } |
84 | ||
108615fc | 85 | AliL3Evaluate::~AliL3Evaluate() |
86 | { | |
c0217812 | 87 | if(fDigitsTree) fDigitsTree->Delete(); |
e1a32fa5 | 88 | if(fDigitsFile) { |
89 | fDigitsFile->Close(); | |
90 | delete fDigitsFile; | |
c0217812 | 91 | } |
c0217812 | 92 | if(fTracks) delete fTracks; |
93 | if(fPtRes) delete fPtRes; | |
94 | if(fNGoodTracksPt) delete fNGoodTracksPt; | |
95 | if(fNFoundTracksPt) delete fNFoundTracksPt; | |
96 | if(fNFakeTracksPt) delete fNFakeTracksPt; | |
97 | if(fTrackEffPt) delete fTrackEffPt; | |
98 | if(fFakeTrackEffPt) delete fFakeTrackEffPt; | |
99 | if(fNGoodTracksEta) delete fNGoodTracksEta; | |
100 | if(fNFoundTracksEta) delete fNFoundTracksEta; | |
101 | if(fNFakeTracksEta) delete fNFakeTracksEta; | |
102 | if(fTrackEffEta) delete fTrackEffEta; | |
103 | if(fFakeTrackEffEta) delete fFakeTrackEffEta; | |
104 | if(fMcIndex) delete [] fMcIndex; | |
105 | if(fMcId) delete [] fMcId; | |
1f79afc0 | 106 | if(fNtuppel) delete fNtuppel; |
108615fc | 107 | } |
108 | ||
bc2f4f0e | 109 | void AliL3Evaluate::Setup(Char_t *trackfile,Char_t *path) |
108615fc | 110 | { |
108615fc | 111 | //Read in the hit and track information from produced files. |
eeddc64d | 112 | AliL3Transform::Init(path); |
108615fc | 113 | |
23908b9b | 114 | fGoodFound = 0; |
115 | fGoodGen = 0; | |
108615fc | 116 | Char_t fname[256]; |
40896357 | 117 | AliL3FileHandler *clusterfile[36][6]; |
108615fc | 118 | for(Int_t s=fMinSlice; s<=fMaxSlice; s++) |
119 | { | |
40896357 | 120 | for(Int_t p=0; p<6; p++) |
108615fc | 121 | { |
1f79afc0 | 122 | fClusters[s][p] = 0; |
108615fc | 123 | clusterfile[s][p] = new AliL3FileHandler(); |
bc2f4f0e | 124 | sprintf(fname,"%s/points_%d_%d.raw",path,s,p); |
108615fc | 125 | if(!clusterfile[s][p]->SetBinaryInput(fname)) |
126 | { | |
127 | LOG(AliL3Log::kError,"AliL3Evaluation::Setup","File Open") | |
128 | <<"Inputfile "<<fname<<" does not exist"<<ENDLOG; | |
bc2f4f0e | 129 | delete clusterfile[s][p]; |
130 | clusterfile[s][p] = 0; | |
131 | continue; | |
108615fc | 132 | } |
133 | fClusters[s][p] = (AliL3SpacePointData*)clusterfile[s][p]->Allocate(); | |
134 | clusterfile[s][p]->Binary2Memory(fNcl[s][p],fClusters[s][p]); | |
135 | clusterfile[s][p]->CloseBinaryInput(); | |
136 | } | |
137 | } | |
78127e35 | 138 | |
23908b9b | 139 | |
108615fc | 140 | AliL3FileHandler *tfile = new AliL3FileHandler(); |
bc2f4f0e | 141 | if(!tfile->SetBinaryInput(trackfile)){ |
142 | LOG(AliL3Log::kError,"AliL3Evaluation::Setup","File Open") | |
143 | <<"Inputfile "<<trackfile<<" does not exist"<<ENDLOG; | |
144 | return; | |
145 | } | |
108615fc | 146 | fTracks = new AliL3TrackArray(); |
147 | tfile->Binary2TrackArray(fTracks); | |
148 | tfile->CloseBinaryInput(); | |
149 | delete tfile; | |
bc2f4f0e | 150 | } |
151 | ||
40896357 | 152 | void AliL3Evaluate::SetupSlow(Char_t *trackfile,Char_t *digitsfile,Char_t *path) |
bc2f4f0e | 153 | { |
154 | //Setup for using the slow simulator. | |
eeddc64d | 155 | AliL3Transform::Init(path); |
bc2f4f0e | 156 | |
40896357 | 157 | fDigitsFile = new TFile(digitsfile,"READ"); |
158 | if(!fDigitsFile->IsOpen()) | |
159 | { | |
160 | LOG(AliL3Log::kError,"AliL3Evaluate::SetupSlow","File Open") | |
161 | <<"Inputfile "<<digitsfile<<" does not exist"<<ENDLOG; | |
162 | return; | |
163 | } | |
bc2f4f0e | 164 | fIsSlow = true; |
165 | Setup(trackfile,path); | |
108615fc | 166 | |
e1a32fa5 | 167 | if(!InitMC()) |
bc2f4f0e | 168 | LOG(AliL3Log::kError,"AliL3Evaluation::SetupSlow","Digits Tree") |
169 | <<"Error setting up digits tree"<<ENDLOG; | |
e1a32fa5 | 170 | |
108615fc | 171 | } |
172 | ||
bc2f4f0e | 173 | void AliL3Evaluate::SetupFast(Char_t *trackfile,Char_t *mcClusterfile,Char_t *path) |
108615fc | 174 | { |
175 | //Setup for using the fast simulator. | |
eeddc64d | 176 | AliL3Transform::Init(path); |
108615fc | 177 | |
178 | fIsSlow = false; | |
c0217812 | 179 | GetFastClusterIDs(path); |
108615fc | 180 | |
181 | fMCclusterfile = new TFile(mcClusterfile); | |
182 | if(!fMCclusterfile->IsOpen()) | |
183 | LOG(AliL3Log::kError,"AliL3Evaluation::SetupFast","File Open") | |
184 | <<"Inputfile "<<mcClusterfile<<" does not exist"<<ENDLOG; | |
185 | ||
bc2f4f0e | 186 | Setup(trackfile,path); |
e1a32fa5 | 187 | InitMC(); |
108615fc | 188 | } |
189 | ||
e1a32fa5 | 190 | Bool_t AliL3Evaluate::InitMC() |
108615fc | 191 | { |
e1a32fa5 | 192 | if(fIsSlow) |
108615fc | 193 | { |
5e0f9911 | 194 | Char_t dname[100]; |
195 | sprintf(dname,"TreeD_%s_0",AliL3Transform::GetParamName()); | |
e1a32fa5 | 196 | fDigitsFile->cd(); |
5e0f9911 | 197 | fDigitsTree = (TTree*)fDigitsFile->Get(dname); |
40896357 | 198 | if(!fDigitsTree) |
199 | { | |
200 | LOG(AliL3Log::kError,"AliL3Evaluate::InitMC","Digits Tree") | |
201 | <<AliL3Log::kHex<<"Error getting digitstree "<<(Int_t)fDigitsTree<<ENDLOG; | |
202 | return false; | |
203 | } | |
e1a32fa5 | 204 | fDigitsTree->GetBranch("Segment")->SetAddress(&fDigits); |
205 | for(Int_t i=0; i<fDigitsTree->GetEntries(); i++) | |
206 | { | |
207 | if(!fDigitsTree->GetEvent(i)) continue; | |
208 | Int_t se,ro,slice,slicerow; | |
209 | fParam->AdjustSectorRow(fDigits->GetID(),se,ro); | |
494fad94 | 210 | AliL3Transform::Sector2Slice(slice,slicerow,se,ro); |
e1a32fa5 | 211 | fRowid[slice][slicerow] = i; |
212 | } | |
108615fc | 213 | } |
214 | ||
e1a32fa5 | 215 | fEventFile->cd(); |
40896357 | 216 | gAlice = (AliRun*)fEventFile->Get("gAlice"); |
108615fc | 217 | if (!gAlice) |
218 | { | |
40896357 | 219 | LOG(AliL3Log::kError,"AliL3Evaluate::InitMC","gAlice") |
108615fc | 220 | <<"AliRun object non existing on file"<<ENDLOG; |
221 | return false; | |
222 | } | |
223 | ||
224 | gAlice->GetEvent(0); | |
e1a32fa5 | 225 | |
108615fc | 226 | return true; |
e1a32fa5 | 227 | |
108615fc | 228 | } |
229 | ||
40896357 | 230 | void AliL3Evaluate::DefineGoodTracks(Int_t *slice,Int_t *padrow,Int_t good_number,Char_t *fname) |
108615fc | 231 | { |
232 | //Loop over MC particles, and mark the good ones | |
233 | //(which the tracker should find...) | |
e1a32fa5 | 234 | |
40896357 | 235 | |
236 | ifstream in(fname); | |
237 | if(in) | |
238 | { | |
239 | LOG(AliL3Log::kInformational,"AliL3Evaluate::DefineGoodTracks","infile") | |
240 | <<"Reading good particles from file "<<fname<<ENDLOG; | |
241 | while (in>>fGoodTracks[fGoodGen].label>>fGoodTracks[fGoodGen].code>> | |
242 | fGoodTracks[fGoodGen].px>>fGoodTracks[fGoodGen].py>>fGoodTracks[fGoodGen].pz>> | |
243 | fGoodTracks[fGoodGen].pt) | |
244 | { | |
245 | fGoodGen++; | |
246 | if (fGoodGen==15000) | |
247 | { | |
248 | LOG(AliL3Log::kWarning,"AliL3Evaluate::DefineGoodTracks","fGoodGen") | |
249 | <<"Too many good tracks"<<ENDLOG; | |
250 | break; | |
251 | } | |
252 | } | |
253 | if (!in.eof()) cerr<<"Read error (good_tracks_tpc) !\n"; | |
254 | return; | |
255 | } | |
256 | ||
257 | ||
108615fc | 258 | AliTPC *TPC = (AliTPC*)gAlice->GetDetector("TPC"); |
e1a32fa5 | 259 | |
108615fc | 260 | TPC->SetParam(fParam); |
e1a32fa5 | 261 | |
108615fc | 262 | Int_t ver = TPC->IsVersion(); |
263 | LOG(AliL3Log::kInformational,"AliL3Evaluate::DefineGoodTracks","TPC version") | |
264 | <<"TPC version "<<ver<<" found on file"<<ENDLOG; | |
265 | ||
108615fc | 266 | Int_t zero=TPC->GetParam()->GetZeroSup(); |
267 | ||
e1a32fa5 | 268 | Int_t np = gAlice->GetNtrack(); |
108615fc | 269 | Int_t *good = new Int_t[np]; |
270 | for(Int_t ii=0; ii<np; ii++) | |
271 | good[ii] = 0; | |
272 | ||
273 | if(ver==1) | |
274 | { | |
275 | if(fIsSlow) | |
276 | LOG(AliL3Log::kError,"AliL3Evaluate::DefineGoodTracks","TPC version") | |
277 | <<"TPC version "<<ver<<" does not match."<<ENDLOG; | |
278 | fMCclusterfile->cd(); | |
279 | AliTPCClustersArray carray; | |
280 | carray.Setup(fParam); | |
281 | carray.SetClusterType("AliTPCcluster"); | |
40896357 | 282 | Char_t cname[100]; |
283 | Int_t eventn = 0; | |
284 | sprintf(cname,"TreeC_TPC_%d",eventn); | |
285 | Bool_t clusterok = carray.ConnectTree(cname); | |
108615fc | 286 | if(!clusterok) |
287 | LOG(AliL3Log::kError,"AliL3Evaluate::DefineGoodTracks","Cluster Array") | |
288 | <<"Error loading clusters from rootfile"<<ENDLOG; | |
78127e35 | 289 | |
bc2f4f0e | 290 | for(Int_t i=0; i<carray.GetTree()->GetEntries(); i++) |
108615fc | 291 | { |
292 | Int_t sec,row,sl,lr; | |
293 | AliSegmentID *s = carray.LoadEntry(i); | |
294 | fParam->AdjustSectorRow(s->GetID(),sec,row); | |
494fad94 | 295 | AliL3Transform::Sector2Slice(sl,lr,sec,row); |
108615fc | 296 | |
e1a32fa5 | 297 | if(sl != slice[0]) {carray.ClearRow(sec,row); continue;} |
bc2f4f0e | 298 | if(lr < padrow[0]) {carray.ClearRow(sec,row); continue;} |
299 | if(lr > padrow[1]) {carray.ClearRow(sec,row); continue;} | |
108615fc | 300 | AliTPCClustersRow *cRow = carray.GetRow(sec,row); |
301 | for(Int_t j=0; j<cRow->GetArray()->GetEntriesFast(); j++) | |
302 | { | |
303 | AliTPCcluster *cluster=(AliTPCcluster*)(*cRow)[j]; | |
304 | Int_t lab=cluster->GetLabel(0); | |
305 | if(lab<0) continue; | |
306 | lab=TMath::Abs(lab); | |
307 | good[lab]++; | |
308 | } | |
309 | if(carray.GetRow(sec,row)) | |
310 | carray.ClearRow(sec,row); | |
311 | } | |
312 | } | |
313 | else if(ver==2) | |
314 | { | |
315 | if(!fIsSlow) | |
316 | LOG(AliL3Log::kError,"AliL3Evaluate::DefineGoodTracks","TPC version") | |
317 | <<"TPC version "<<ver<<" does not match."<<ENDLOG; | |
318 | Int_t *count = new Int_t[np]; //np number of particles. | |
319 | Int_t i; | |
355debe1 | 320 | Float_t xyz[3]; |
108615fc | 321 | for (i=0; i<np; i++) count[i]=0; |
e1a32fa5 | 322 | for(Int_t sl=slice[0]; sl<=slice[1]; sl++) |
323 | { | |
324 | for (i=padrow[0]; i<=padrow[1]; i++) { | |
325 | Int_t index = fRowid[sl][i]; | |
326 | if (!fDigitsTree->GetEvent(index)) continue; | |
327 | Int_t sec,row; | |
328 | fParam->AdjustSectorRow(fDigits->GetID(),sec,row); | |
329 | fDigits->First(); | |
330 | do { | |
331 | Int_t it=fDigits->CurrentRow(), ip=fDigits->CurrentColumn(); | |
332 | Short_t dig = fDigits->GetDigit(it,ip); | |
355debe1 | 333 | |
334 | if(dig<=fParam->GetZeroSup()) continue; | |
e8541476 | 335 | /* |
336 | if(it < fParam->GetMaxTBin()-1 && it > 0) | |
355debe1 | 337 | if(fDigits->GetDigit(it+1,ip) <= fParam->GetZeroSup() |
e8541476 | 338 | && fDigits->GetDigit(it-1,ip) <= fParam->GetZeroSup()) |
339 | continue; | |
340 | */ | |
494fad94 | 341 | AliL3Transform::Raw2Local(xyz,sec,row,ip,it); |
355debe1 | 342 | if(fParam->GetPadRowRadii(sec,row)<230./250.*fabs(xyz[2])) |
343 | continue; | |
344 | ||
e1a32fa5 | 345 | Int_t idx0=fDigits->GetTrackID(it,ip,0); |
346 | Int_t idx1=fDigits->GetTrackID(it,ip,1); | |
347 | Int_t idx2=fDigits->GetTrackID(it,ip,2); | |
e170e0f9 | 348 | |
e1a32fa5 | 349 | if (idx0>=0 && dig>=zero) count[idx0]+=1; |
350 | if (idx1>=0 && dig>=zero) count[idx1]+=1; | |
351 | if (idx2>=0 && dig>=zero) count[idx2]+=1; | |
352 | } while (fDigits->Next()); | |
353 | ||
354 | for (Int_t j=0; j<np; j++) | |
355 | { | |
356 | if (count[j]>1) //at least two digits at this padrow | |
357 | good[j]++; | |
358 | ||
359 | count[j]=0; | |
360 | } | |
108615fc | 361 | } |
108615fc | 362 | } |
108615fc | 363 | delete[] count; |
364 | } | |
e1a32fa5 | 365 | |
108615fc | 366 | else |
367 | { | |
368 | LOG(AliL3Log::kError,"AliL3Evaluation::FillEffHistos","TPC version") | |
369 | <<"No valid TPC version found"<<ENDLOG; | |
e1a32fa5 | 370 | return; |
108615fc | 371 | } |
372 | ||
e1a32fa5 | 373 | |
374 | //The following code has been taken from offlinemacro->AliTPCComparison.C | |
375 | ||
376 | TTree *TH=gAlice->TreeH(); | |
377 | Int_t npart=(Int_t)TH->GetEntries(); | |
e1a32fa5 | 378 | Int_t max = 15000; |
379 | while (npart--) { | |
380 | AliTPChit *hit0=0; | |
381 | ||
40896357 | 382 | TPC->ResetHits(); |
383 | TH->GetEvent(npart); | |
384 | AliTPChit * hit = (AliTPChit*) TPC->FirstHit(-1); | |
385 | while (hit){ | |
386 | if (hit->fQ==0.) break; | |
387 | hit = (AliTPChit*) TPC->NextHit(); | |
388 | } | |
389 | if (hit) { | |
390 | hit0 = new AliTPChit(*hit); //Make copy of hit | |
391 | hit = hit0; | |
392 | } | |
393 | else continue; | |
394 | AliTPChit *hit1=(AliTPChit*)TPC->NextHit(); | |
395 | if (hit1==0) continue; | |
396 | if (hit1->fQ != 0.) continue; | |
397 | Int_t i=hit->Track(); | |
398 | TParticle *p = (TParticle*)gAlice->Particle(i); | |
399 | ||
400 | //printf("Checking particle %d with code %d\n",i,p->GetPdgCode()); | |
401 | if (p->GetFirstMother()>=0) continue; //secondary particle | |
355debe1 | 402 | if (good[i] < good_number || good[i]>176) continue; |
40896357 | 403 | if (p->Pt()<fMinGoodPt) continue; |
404 | if (TMath::Abs(p->Pz()/p->Pt())>0.999) continue; | |
405 | printf("Found good particle %d, nHits %d\n",i,good[i]); | |
406 | ||
407 | fGoodTracks[fGoodGen].label=i; | |
408 | fGoodTracks[fGoodGen].code=p->GetPdgCode(); | |
409 | //**** px py pz - in global coordinate system, x y z - in local ! | |
410 | fGoodTracks[fGoodGen].px=hit->X(); fGoodTracks[fGoodGen].py=hit->Y(); fGoodTracks[fGoodGen].pz=hit->Z(); | |
411 | fGoodTracks[fGoodGen].pt = p->Pt(); | |
412 | fGoodGen++; | |
413 | //nt++; | |
414 | if (hit0) delete hit0; | |
415 | if (fGoodGen==max) {cerr<<"Too many good tracks !n"; break;} | |
416 | } | |
e170e0f9 | 417 | |
108615fc | 418 | delete [] good; |
40896357 | 419 | |
420 | LOG(AliL3Log::kInformational,"AliL3Evaluate::DefineGoodParticles","Eval") | |
421 | <<AliL3Log::kDec<<"Found "<<fGoodGen<<"good tracks"<<ENDLOG; | |
422 | ||
423 | ofstream out(fname); | |
424 | if(out) | |
425 | { | |
426 | for (Int_t ngd=0; ngd<fGoodGen; ngd++) | |
427 | out<<fGoodTracks[ngd].label<<' '<<fGoodTracks[ngd].code<<' '<< | |
428 | fGoodTracks[ngd].px<<' '<<fGoodTracks[ngd].py<<' '<<fGoodTracks[ngd].pz<<' '<< | |
429 | fGoodTracks[ngd].pt<<endl; | |
430 | ||
431 | } | |
432 | else | |
433 | LOG(AliL3Log::kError,"AliL3Evaluate::DefineGoodParticles","file") | |
434 | <<"Can not open file containing good tracks"<<ENDLOG; | |
435 | ||
436 | out.close(); | |
437 | ||
108615fc | 438 | } |
439 | ||
440 | void AliL3Evaluate::EvaluatePatch(Int_t slice,Int_t patch,Int_t min_points,Int_t good_number) | |
441 | { | |
442 | //Make efficiency plots for tracking on patch level (before any merging). | |
443 | ||
40896357 | 444 | Int_t row[6][2] = {{ 0, 45},{46,77},{78,109},{110,141},{142,175}}; |
e1a32fa5 | 445 | Int_t sl[2] ={slice,slice}; |
446 | DefineGoodTracks(sl,row[patch],good_number); | |
108615fc | 447 | SetMinPoints(min_points); |
448 | AssignIDs(); | |
bc2f4f0e | 449 | CreateHistos(); |
e1a32fa5 | 450 | FillEffHistos(); |
bc2f4f0e | 451 | CalcEffHistos(); |
108615fc | 452 | } |
453 | ||
454 | void AliL3Evaluate::EvaluateSlice(Int_t slice,Int_t min_points,Int_t good_number) | |
455 | { | |
456 | //Make efficiency plots for tracking on a slice (after merging). | |
457 | //min_points = minimum points on track to be considered for evaluation | |
458 | //good_number = minimum hits (padrows) produced by simulated track for consideration. | |
459 | ||
ae97a0b9 | 460 | Int_t row[2] = {0,175}; |
e1a32fa5 | 461 | Int_t sl[2] ={slice,slice}; |
462 | DefineGoodTracks(sl,row,good_number); | |
23908b9b | 463 | |
108615fc | 464 | SetMinPoints(min_points); |
23908b9b | 465 | |
108615fc | 466 | AssignIDs(); |
bc2f4f0e | 467 | CreateHistos(); |
e1a32fa5 | 468 | FillEffHistos(); |
bc2f4f0e | 469 | CalcEffHistos(); |
108615fc | 470 | } |
471 | ||
40896357 | 472 | void AliL3Evaluate::EvaluateGlobal(Int_t min_points,Int_t good_number,Char_t *fname) |
108615fc | 473 | { |
474 | //Make efficiency plots for tracking on several slices. | |
bc2f4f0e | 475 | |
e8541476 | 476 | Int_t row[2] = {AliL3Transform::GetFirstRow(0),AliL3Transform::GetLastRow(5)}; |
bc2f4f0e | 477 | SetMinPoints(min_points); |
e1a32fa5 | 478 | Int_t slice[2] = {fMinSlice,fMaxSlice}; |
40896357 | 479 | DefineGoodTracks(slice,row,good_number,fname); |
bc2f4f0e | 480 | AssignIDs(); |
355debe1 | 481 | CreateHistos(5,0,2); |
40896357 | 482 | printf("filling histos\n"); |
e1a32fa5 | 483 | FillEffHistos(); |
40896357 | 484 | printf("done fillling\n"); |
bc2f4f0e | 485 | CalcEffHistos(); |
108615fc | 486 | } |
487 | ||
108615fc | 488 | void AliL3Evaluate::AssignIDs() |
489 | { | |
490 | //Assign MC id to the tracks. | |
491 | ||
108615fc | 492 | fTracks->QSort(); |
493 | LOG(AliL3Log::kDebug,"AliL3Evaluate::AssignIDs","Track Loop") | |
494 | <<"Assigning MC id to the found tracks...."<<ENDLOG; | |
495 | for(Int_t i=0; i<fTracks->GetNTracks(); i++) | |
496 | { | |
497 | AliL3Track *track = (AliL3Track*)fTracks->GetCheckedTrack(i); | |
c0217812 | 498 | if(!track) continue; |
108615fc | 499 | if(track->GetNumberOfPoints() < fMinPointsOnTrack) break; |
23908b9b | 500 | |
23908b9b | 501 | fGoodFound++; |
c0217812 | 502 | Int_t tID = GetMCTrackLabel(track); |
108615fc | 503 | track->SetMCid(tID); |
e170e0f9 | 504 | printf("track %i id %d nHits %d\n",i,tID,track->GetNumberOfPoints()); |
108615fc | 505 | } |
108615fc | 506 | } |
507 | ||
508 | ||
509 | struct S {Int_t lab; Int_t max;}; | |
c0217812 | 510 | Int_t AliL3Evaluate::GetMCTrackLabel(AliL3Track *track){ |
108615fc | 511 | //Returns the MCtrackID of the belonging clusters. |
512 | //If MCLabel < 0, means that track is fake. | |
513 | //Fake track means that more than 10 percent of clusters are assigned incorrectly. | |
494fad94 | 514 | |
515 | #ifdef do_mc | |
108615fc | 516 | Int_t num_of_clusters = track->GetNumberOfPoints(); |
517 | S *s=new S[num_of_clusters]; | |
518 | Int_t i; | |
519 | for (i=0; i<num_of_clusters; i++) s[i].lab=s[i].max=0; | |
520 | UInt_t *hitnum = track->GetHitNumbers(); | |
521 | UInt_t id; | |
108615fc | 522 | |
355debe1 | 523 | //Int_t **trackID = GetClusterIDs(track); |
108615fc | 524 | |
525 | Int_t lab=123456789; | |
526 | for (i=0; i<num_of_clusters; i++) | |
527 | { | |
528 | //Tricks to get the clusters belonging to this track: | |
529 | id = hitnum[i]; | |
530 | Int_t slice = (id>>25) & 0x7f; | |
531 | Int_t patch = (id>>22) & 0x7; | |
532 | UInt_t pos = id&0x3fffff; | |
533 | ||
534 | AliL3SpacePointData *points = fClusters[slice][patch]; | |
c0217812 | 535 | if(!points) continue; |
108615fc | 536 | if(pos>=fNcl[slice][patch]) |
537 | { | |
538 | LOG(AliL3Log::kError,"AliL3Evaluate::GetMCTrackLabel","Clusterarray") | |
539 | <<AliL3Log::kDec<<"ERROR"<<ENDLOG; | |
540 | continue; | |
541 | } | |
542 | ||
543 | //Get the label of the cluster: | |
355debe1 | 544 | |
545 | lab=abs(points[pos].fTrackID[0]); | |
108615fc | 546 | Int_t j; |
547 | for (j=0; j<num_of_clusters; j++) | |
548 | if (s[j].lab==lab || s[j].max==0) break; | |
549 | s[j].lab=lab; | |
550 | s[j].max++; | |
551 | } | |
552 | ||
553 | Int_t max=0; | |
554 | for (i=0; i<num_of_clusters; i++) | |
555 | if (s[i].max>max) {max=s[i].max; lab=s[i].lab;} | |
556 | ||
557 | delete[] s; | |
558 | ||
559 | for (i=0; i<num_of_clusters; i++) | |
560 | { | |
561 | id = hitnum[i]; | |
562 | Int_t slice = (id>>25) & 0x7f; | |
563 | Int_t patch = (id>>22) & 0x7; | |
564 | UInt_t pos = id&0x3fffff; | |
565 | ||
566 | AliL3SpacePointData *points = fClusters[slice][patch]; | |
c0217812 | 567 | if(!points) continue; |
108615fc | 568 | if(pos>=fNcl[slice][patch]) |
569 | { | |
570 | LOG(AliL3Log::kError,"AliL3Evaluate::GetMCTrackLabel","Clusterarray") | |
571 | <<AliL3Log::kDec<<"ERROR"<<ENDLOG; | |
572 | continue; | |
573 | } | |
574 | ||
355debe1 | 575 | if (abs(points[pos].fTrackID[1]) == lab || |
576 | abs(points[pos].fTrackID[2]) == lab ) max++; | |
108615fc | 577 | } |
578 | ||
579 | //check if more than 10% of the clusters are incorrectly assigned (fake track): | |
580 | if (1.-Float_t(max)/num_of_clusters > 0.10) | |
581 | { | |
582 | return -lab; | |
583 | } | |
584 | ||
108615fc | 585 | return lab; |
494fad94 | 586 | #else //If we are running with mc_ids or not |
587 | return 0; | |
588 | #endif | |
589 | ||
108615fc | 590 | } |
591 | ||
c0217812 | 592 | Int_t **AliL3Evaluate::GetClusterIDs(AliL3Track *track) |
108615fc | 593 | { |
594 | //Return the MC information of all clusters belonging to track. | |
c0217812 | 595 | |
108615fc | 596 | Int_t num_of_clusters = track->GetNumberOfPoints(); |
597 | Int_t **trackID = new Int_t*[num_of_clusters]; | |
598 | ||
599 | UInt_t *hitnum = track->GetHitNumbers(); | |
600 | UInt_t id; | |
108615fc | 601 | |
602 | Float_t xyz[3]; | |
c0217812 | 603 | Int_t padrow; |
108615fc | 604 | for(Int_t i=0; i<num_of_clusters; i++) |
605 | { | |
606 | id = hitnum[i]; | |
607 | Int_t slice = (id>>25) & 0x7f; | |
608 | Int_t patch = (id>>22) & 0x7; | |
609 | UInt_t pos = id&0x3fffff; | |
610 | ||
611 | AliL3SpacePointData *points = fClusters[slice][patch]; | |
612 | ||
613 | if(!points) | |
614 | { | |
615 | LOG(AliL3Log::kError,"AliL3Evaluate::GetClusterIDs","Clusterarray") | |
616 | <<"No points at slice "<<slice<<" patch "<<patch<<" pos "<<pos<<ENDLOG; | |
617 | continue; | |
618 | } | |
619 | if(pos>=fNcl[slice][patch]) | |
620 | { | |
621 | LOG(AliL3Log::kError,"AliL3Evaluate::GetClusterIDs","Clusterarray") | |
622 | <<AliL3Log::kDec<<"ERROR"<<ENDLOG; | |
623 | continue; | |
624 | } | |
625 | ||
626 | xyz[0] = points[pos].fX; | |
627 | xyz[1] = points[pos].fY; | |
628 | xyz[2] = points[pos].fZ; | |
629 | //sector = points[pos].fSector; | |
630 | padrow = points[pos].fPadRow; | |
631 | Int_t se,ro; | |
494fad94 | 632 | AliL3Transform::Slice2Sector(slice,padrow,se,ro); |
633 | AliL3Transform::Global2Raw(xyz,se,ro); | |
108615fc | 634 | |
635 | if(fIsSlow) | |
636 | { | |
637 | Int_t p = fRowid[slice][padrow]; | |
108615fc | 638 | if(!fDigitsTree->GetEvent(p)) |
639 | LOG(AliL3Log::kError,"AliL3Evaluate::GetClusterIDs","Digits Tree") | |
640 | <<"Error reading digits tree"<<ENDLOG; | |
641 | ||
642 | trackID[i] = new Int_t[3]; | |
e170e0f9 | 643 | trackID[i][0] = fDigits->GetTrackID((Int_t)rint(xyz[2]),(Int_t)rint(xyz[1]),0); |
644 | trackID[i][1] = fDigits->GetTrackID((Int_t)rint(xyz[2]),(Int_t)rint(xyz[1]),1); | |
645 | trackID[i][2] = fDigits->GetTrackID((Int_t)rint(xyz[2]),(Int_t)rint(xyz[1]),2); | |
646 | //if(trackID[i][0]==6 || trackID[i][0]==32) | |
647 | //printf("trackID %d, padrow %d pad %d time %d\n",fDigits->GetTrackID((Int_t)rint(xyz[2]),(Int_t)rint(xyz[1]),0),padrow,(int)rint(xyz[1]),(int)rint(xyz[2])); | |
648 | /*if(trackID[i][0]<0) | |
649 | { | |
650 | ||
651 | printf("trackID %d, padrow %d pad %d time %d\n",trackID[i][0],padrow,(int)rint(xyz[1]),(int)rint(xyz[2])); | |
652 | printf("on the side %d %d %d %d\n",fDigits->GetTrackID(((int)rint(xyz[2])-1),((int)rint(xyz[1])),0),fDigits->GetTrackID(((int)rint(xyz[2])+1),((int)rint(xyz[1])),0),fDigits->GetTrackID(((int)rint(xyz[2])),((int)rint(xyz[1])-1),0),fDigits->GetTrackID(((int)rint(xyz[2])),((int)rint(xyz[1])+1),0)); | |
653 | }*/ | |
108615fc | 654 | } |
655 | else | |
656 | { | |
657 | Int_t tmp_pid=0; | |
c0217812 | 658 | for(Int_t ii=0; ii<fNFastPoints; ii++) |
108615fc | 659 | { |
c0217812 | 660 | tmp_pid = fMcId[ii]; |
661 | if(fMcIndex[ii] == id) break; | |
108615fc | 662 | } |
663 | trackID[i] = new Int_t[3]; | |
664 | trackID[i][0] = tmp_pid; | |
665 | trackID[i][1] = -1; | |
666 | trackID[i][2] = -1; | |
667 | } | |
668 | } | |
e170e0f9 | 669 | |
108615fc | 670 | return trackID; |
671 | } | |
672 | ||
c0217812 | 673 | void AliL3Evaluate::GetFastClusterIDs(Char_t *path) |
108615fc | 674 | { |
675 | //Get the MC id of space points in case of using the fast simulator. | |
c0217812 | 676 | char fname[256]; |
677 | sprintf(fname,"%s/point_mc.dat",path); | |
678 | FILE *infile = fopen(fname,"r"); | |
679 | if(!infile) return; | |
108615fc | 680 | Int_t hitid,hitmc,i; |
681 | ||
682 | for(i=0; ; i++) | |
683 | if(fscanf(infile,"%d %d",&hitid,&hitmc)==EOF) break; | |
108615fc | 684 | rewind(infile); |
c0217812 | 685 | fNFastPoints = i; |
686 | fMcId = new Int_t[fNFastPoints]; | |
687 | fMcIndex = new UInt_t[fNFastPoints]; | |
108615fc | 688 | |
c0217812 | 689 | for(i=0; i<fNFastPoints; i++) |
108615fc | 690 | { |
691 | if(fscanf(infile,"%d %d",&hitid,&hitmc)==EOF) break; | |
c0217812 | 692 | fMcId[i] = hitmc; |
693 | fMcIndex[i] = hitid; | |
108615fc | 694 | } |
695 | fclose(infile); | |
c0217812 | 696 | } |
697 | ||
698 | void AliL3Evaluate::CreateHistos(Int_t nbin,Int_t xlow,Int_t xup) | |
699 | { | |
700 | //Create the histograms | |
701 | ||
40896357 | 702 | LOG(AliL3Log::kInformational,"AliL3Evaluate::CreateHistos","Allocating") |
703 | <<"Creating histograms..."<<ENDLOG; | |
704 | ||
1f79afc0 | 705 | fNtuppel = new TNtuple("fNtuppel","Pt resolution","pt_gen:pt_found:nHits"); |
706 | ||
c0217812 | 707 | fPtRes = new TH1F("fPtRes","Relative Pt resolution",30,-10.,10.); |
708 | fNGoodTracksPt = new TH1F("fNGoodTracksPt","Good tracks vs pt",nbin,xlow,xup); | |
709 | fNFoundTracksPt = new TH1F("fNFoundTracksPt","Found tracks vs pt",nbin,xlow,xup); | |
710 | fNFakeTracksPt = new TH1F("fNFakeTracksPt","Fake tracks vs pt",nbin,xlow,xup); | |
711 | fTrackEffPt = new TH1F("fTrackEffPt","Tracking efficiency vs pt",nbin,xlow,xup); | |
712 | fFakeTrackEffPt = new TH1F("fFakeTrackEffPt","Efficiency for fake tracks vs pt",nbin,xlow,xup); | |
713 | ||
714 | fNGoodTracksEta = new TH1F("fNGoodTracksEta","Good tracks vs eta",20,-50,50); | |
715 | fNFoundTracksEta = new TH1F("fNFoundTracksEta","Found tracks vs eta",20,-50,50); | |
716 | fNFakeTracksEta = new TH1F("fNFakeTracksEta","Fake tracks vs eta",20,-50,50); | |
717 | fTrackEffEta = new TH1F("fTrackEffEta","Tracking efficienct vs eta",20,-50,50); | |
718 | fFakeTrackEffEta = new TH1F("fFakeTrackEffEta","Efficiency for fake tracks vs eta",20,-50,50); | |
40896357 | 719 | printf("finished creating histos\n"); |
c0217812 | 720 | } |
eeddc64d | 721 | |
355debe1 | 722 | /* |
e1a32fa5 | 723 | void AliL3Evaluate::FillEffHistos() |
c0217812 | 724 | { |
725 | //Fill the efficiency histograms. | |
355debe1 | 726 | |
e1a32fa5 | 727 | for(Int_t i=0; i<fGoodGen; i++) |
c0217812 | 728 | { |
e1a32fa5 | 729 | Double_t ptg=fGoodTracks[i].pt,pzg=fGoodTracks[i].pz; |
c0217812 | 730 | Float_t dipangle=TMath::ATan2(pzg,ptg)*180./TMath::Pi(); |
355debe1 | 731 | //printf("filling particle with pt %f and dipangle %f\n",ptg,dipangle); |
c0217812 | 732 | fNGoodTracksPt->Fill(ptg); |
733 | fNGoodTracksEta->Fill(dipangle); | |
734 | Int_t found = 0; | |
735 | for(Int_t k=0; k<fTracks->GetNTracks(); k++) | |
736 | { | |
737 | AliL3Track *track = fTracks->GetCheckedTrack(k); | |
738 | if(!track) continue; | |
739 | Int_t nHits = track->GetNumberOfPoints(); | |
740 | if(nHits < fMinPointsOnTrack) break; | |
741 | ||
742 | Int_t tracklabel; | |
743 | tracklabel = track->GetMCid(); | |
355debe1 | 744 | //printf("evaluating track id %d\n",tracklabel); |
e1a32fa5 | 745 | if(TMath::Abs(tracklabel) != fGoodTracks[i].label) continue; |
c0217812 | 746 | found=1; |
e1a32fa5 | 747 | if(tracklabel == fGoodTracks[i].label) {fNFoundTracksPt->Fill(ptg); fNFoundTracksEta->Fill(dipangle);} |
c0217812 | 748 | else {fNFakeTracksPt->Fill(ptg); fNFakeTracksEta->Fill(dipangle);} |
749 | Float_t pt=track->GetPt(); | |
750 | fPtRes->Fill((pt-ptg)/ptg*100.); | |
1f79afc0 | 751 | fNtuppel->Fill(ptg,pt,nHits); |
c0217812 | 752 | break; |
753 | ||
754 | } | |
755 | } | |
756 | } | |
355debe1 | 757 | */ |
eeddc64d | 758 | |
355debe1 | 759 | void AliL3Evaluate::FillEffHistos() |
760 | { | |
761 | //Fill the efficiency histograms. | |
762 | ||
763 | cout<<endl<<"Note: Doing NAIVE evaluation "<<endl; | |
764 | for(Int_t i=0; i<fGoodGen; i++) | |
765 | { | |
766 | Double_t ptg=fGoodTracks[i].pt,pzg=fGoodTracks[i].pz; | |
767 | Float_t dipangle=TMath::ATan2(pzg,ptg)*180./TMath::Pi(); | |
768 | //printf("filling particle with pt %f and dipangle %f\n",ptg,dipangle); | |
769 | fNGoodTracksPt->Fill(ptg); | |
770 | fNGoodTracksEta->Fill(dipangle); | |
771 | ||
772 | } | |
773 | ||
774 | for(Int_t k=0; k<fTracks->GetNTracks(); k++) | |
775 | { | |
776 | AliL3Track *track = fTracks->GetCheckedTrack(k); | |
777 | if(!track) continue; | |
778 | Int_t nHits = track->GetNumberOfPoints(); | |
779 | if(nHits < fMinPointsOnTrack) break; | |
780 | if(track->GetPt()<fMinGoodPt) continue; | |
781 | if(fabs(track->GetPseudoRapidity())>0.9) continue; | |
782 | ||
783 | fNFoundTracksPt->Fill(track->GetPt()); fNFoundTracksEta->Fill(track->GetPseudoRapidity()); | |
784 | //Float_t pt=track->GetPt(); | |
785 | //fPtRes->Fill((pt-ptg)/ptg*100.); | |
786 | //fNtuppel->Fill(ptg,pt,nHits); | |
787 | ||
788 | } | |
789 | } | |
c0217812 | 790 | |
791 | void AliL3Evaluate::CalcEffHistos(){ | |
792 | ||
793 | Stat_t ngood=fNGoodTracksPt->GetEntries(); | |
794 | Stat_t nfound=fNFoundTracksPt->GetEntries(); | |
795 | Stat_t nfake=fNFakeTracksPt->GetEntries(); | |
796 | ||
797 | LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency") | |
798 | <<AliL3Log::kDec<<"There was "<<ngood<<" generated good tracks"<<ENDLOG; | |
799 | LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency") | |
800 | <<AliL3Log::kDec<<"Found "<<nfound<<" tracks"<<ENDLOG; | |
801 | LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency") | |
802 | <<AliL3Log::kDec<<"Integral efficiency is about "<<nfound/ngood*100<<ENDLOG; | |
803 | LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency") | |
804 | <<AliL3Log::kDec<<"Fake tracks relative is about "<<nfake/ngood*100<<ENDLOG; | |
805 | ||
355debe1 | 806 | LOG(AliL3Log::kInformational,"AliL3Evaluate::FillEffHistos","Efficiency") |
807 | <<AliL3Log::kDec<<"Naive efficiency "<<(Double_t)fGoodFound/(Double_t)fGoodGen<<ENDLOG; | |
808 | ||
809 | ||
c0217812 | 810 | fNFoundTracksPt->Sumw2(); fNGoodTracksPt->Sumw2(); |
811 | fTrackEffPt->Divide(fNFoundTracksPt,fNGoodTracksPt,1,1.,"b"); | |
812 | fFakeTrackEffPt->Divide(fNFakeTracksPt,fNGoodTracksPt,1,1.,"b"); | |
813 | fTrackEffPt->SetMaximum(1.4); | |
814 | fTrackEffPt->SetXTitle("P_{T} [GeV]"); | |
815 | fTrackEffPt->SetLineWidth(2); | |
816 | fFakeTrackEffPt->SetFillStyle(3013); | |
817 | fTrackEffPt->SetLineColor(4); | |
818 | fFakeTrackEffPt->SetFillColor(2); | |
819 | ||
820 | fNFoundTracksEta->Sumw2(); fNGoodTracksEta->Sumw2(); | |
821 | fTrackEffEta->Divide(fNFoundTracksEta,fNGoodTracksEta,1,1.,"b"); | |
822 | fFakeTrackEffEta->Divide(fNFakeTracksEta,fNGoodTracksEta,1,1.,"b"); | |
823 | fTrackEffEta->SetMaximum(1.4); | |
824 | fTrackEffEta->SetXTitle("#lambda [degrees]"); | |
825 | fTrackEffEta->SetLineWidth(2); | |
826 | fFakeTrackEffEta->SetFillStyle(3013); | |
827 | fTrackEffEta->SetLineColor(4); | |
828 | fFakeTrackEffEta->SetFillColor(2); | |
829 | ||
108615fc | 830 | } |
831 | ||
832 | void AliL3Evaluate::Write2File(Char_t *outputfile) | |
833 | { | |
834 | //Write histograms to file: | |
835 | ||
836 | TFile *of = new TFile(outputfile,"RECREATE"); | |
837 | if(!of->IsOpen()) | |
838 | { | |
839 | LOG(AliL3Log::kError,"AliL3Evaluate::Write2File","File Open") | |
840 | <<"Problems opening rootfile"<<ENDLOG; | |
841 | return; | |
842 | } | |
843 | ||
844 | of->cd(); | |
1f79afc0 | 845 | fNtuppel->Write(); |
108615fc | 846 | fPtRes->Write(); |
847 | fNGoodTracksPt->Write(); | |
848 | fNFoundTracksPt->Write(); | |
849 | fNFakeTracksPt->Write(); | |
850 | fTrackEffPt->Write(); | |
851 | fFakeTrackEffPt->Write(); | |
852 | fNGoodTracksEta->Write(); | |
853 | fNFoundTracksEta->Write(); | |
854 | fNFakeTracksEta->Write(); | |
855 | fTrackEffEta->Write(); | |
856 | fFakeTrackEffEta->Write(); | |
857 | ||
858 | of->Close(); | |
859 | delete of; | |
860 | ||
861 | } | |
862 | ||
1f79afc0 | 863 | TNtuple *AliL3Evaluate::CalculateResiduals() |
864 | { | |
865 | ||
866 | TNtuple *ntuppel=new TNtuple("ntuppel","Residuals","residual_trans:residual_long:zHit:pt:dipangle:beta:padrow:nHits"); | |
78127e35 | 867 | |
1f79afc0 | 868 | for(int f=fMinSlice; f<=fMaxSlice; f++) |
869 | { | |
870 | AliL3FileHandler *tfile = new AliL3FileHandler(); | |
871 | char fname[256]; | |
872 | sprintf(fname,"tracks_tr_%d_0.raw",f); | |
873 | if(!tfile->SetBinaryInput(fname)){ | |
874 | LOG(AliL3Log::kError,"AliL3Evaluation::Setup","File Open") | |
875 | <<"Inputfile "<<fname<<" does not exist"<<ENDLOG; | |
876 | return 0; | |
877 | } | |
878 | fTracks = new AliL3TrackArray(); | |
879 | tfile->Binary2TrackArray(fTracks); | |
880 | tfile->CloseBinaryInput(); | |
881 | delete tfile; | |
78127e35 | 882 | printf("Looking in slice %d\n",f); |
1f79afc0 | 883 | for(Int_t i=0; i<fTracks->GetNTracks(); i++) |
884 | { | |
885 | ||
886 | AliL3Track *track = (AliL3Track*)fTracks->GetCheckedTrack(i); | |
887 | if(!track) continue; | |
888 | ||
889 | track->CalculateHelix(); | |
890 | UInt_t *hitnum = track->GetHitNumbers(); | |
891 | UInt_t id; | |
892 | ||
893 | Float_t xyz[3]; | |
894 | Int_t padrow; | |
78127e35 | 895 | for(Int_t j=0; j<track->GetNumberOfPoints()-1; j++) |
1f79afc0 | 896 | { |
897 | id = hitnum[j]; | |
898 | Int_t slice = (id>>25) & 0x7f; | |
899 | Int_t patch = (id>>22) & 0x7; | |
900 | UInt_t pos = id&0x3fffff; | |
901 | ||
78127e35 | 902 | //if(slice!=1) continue; |
903 | ||
1f79afc0 | 904 | AliL3SpacePointData *points = fClusters[slice][patch]; |
905 | ||
906 | if(!points) | |
907 | { | |
908 | LOG(AliL3Log::kError,"AliL3Evaluate::CalculateResiduals","Clusterarray") | |
909 | <<"No points at slice "<<slice<<" patch "<<patch<<" pos "<<pos<<ENDLOG; | |
910 | continue; | |
911 | } | |
912 | if(pos>=fNcl[slice][patch]) | |
913 | { | |
914 | LOG(AliL3Log::kError,"AliL3Evaluate::CalculateResiduals","Clusterarray") | |
915 | <<AliL3Log::kDec<<"ERROR"<<ENDLOG; | |
916 | continue; | |
917 | } | |
918 | ||
919 | xyz[0] = points[pos].fX; | |
920 | xyz[1] = points[pos].fY; | |
921 | xyz[2] = points[pos].fZ; | |
922 | padrow = points[pos].fPadRow; | |
494fad94 | 923 | //AliL3Transform::Global2Local(xyz,slice); |
1f79afc0 | 924 | |
925 | Float_t xyz_cross[3]; | |
926 | track->GetCrossingPoint(padrow,xyz_cross); | |
e1a32fa5 | 927 | Double_t beta = track->GetCrossingAngle(padrow); |
1f79afc0 | 928 | |
929 | Double_t yres = xyz_cross[1] - xyz[1]; | |
930 | Double_t zres = xyz_cross[2] - xyz[2]; | |
931 | ||
932 | Double_t dipangle = atan(track->GetTgl()); | |
933 | ntuppel->Fill(yres,zres,xyz_cross[2],track->GetPt(),dipangle,beta,padrow,track->GetNumberOfPoints()); | |
934 | ||
935 | } | |
936 | } | |
937 | if(fTracks) | |
938 | delete fTracks; | |
939 | } | |
940 | return ntuppel; | |
941 | } | |
942 | ||
e1a32fa5 | 943 | TNtuple *AliL3Evaluate::EvaluatePoints(Char_t *rootfile) |
78127e35 | 944 | { |
e1a32fa5 | 945 | //Compare points to the exact crossing points of track and padrows. |
1b76fd18 | 946 | //The input file to this function, contains the exact clusters calculated |
947 | //in AliTPC::Hits2ExactClusters. | |
e1a32fa5 | 948 | |
355debe1 | 949 | cout<<"Evaluating points"<<endl; |
e1a32fa5 | 950 | TNtuple *ntuppel = new TNtuple("ntuppel","residuals","slice:padrow:resy:resz:zHit:pt"); |
78127e35 | 951 | |
355debe1 | 952 | TFile *exfile = TFile::Open(rootfile); |
953 | if(!exfile) | |
954 | { | |
955 | cerr<<"Error opening rootfile "<<rootfile<<endl; | |
956 | return 0; | |
957 | } | |
958 | ||
5e0f9911 | 959 | AliTPCParam *param = (AliTPCParam*)exfile->Get(AliL3Transform::GetParamName()); |
e1a32fa5 | 960 | |
961 | //Get the exact clusters from file: | |
962 | AliTPCClustersArray *arr = new AliTPCClustersArray; | |
963 | arr->Setup(param); | |
964 | arr->SetClusterType("AliComplexCluster"); | |
965 | char treeName[500]; | |
966 | sprintf(treeName,"TreeCExact_%s",param->GetTitle()); | |
294ee877 | 967 | Bool_t clusterok = arr->ConnectTree(treeName);//Segment Tree (for offline clusters) |
e1a32fa5 | 968 | if(!clusterok) {printf("AliL3Evaluate::EvaluatePoints : Error in clusterloading\n"); return 0;} |
969 | ||
355debe1 | 970 | cout<<"Entering loop with "<<(Int_t)arr->GetTree()->GetEntries()<<endl; |
e1a32fa5 | 971 | for(Int_t i=0; i<arr->GetTree()->GetEntries(); i++) |
78127e35 | 972 | { |
e1a32fa5 | 973 | //Get the exact clusters for this row: |
974 | Int_t cursec,currow; | |
975 | AliSegmentID *s = arr->LoadEntry(i); | |
976 | param->AdjustSectorRow(s->GetID(),cursec,currow); | |
78127e35 | 977 | |
e1a32fa5 | 978 | AliTPCClustersRow *ro = (AliTPCClustersRow *)arr->GetRow(cursec,currow); |
979 | TClonesArray *clusters = ro->GetArray(); | |
980 | int num_of_offline=clusters->GetEntriesFast(); | |
78127e35 | 981 | |
e1a32fa5 | 982 | //Get the found clusters: |
983 | Int_t slice,padrow; | |
494fad94 | 984 | AliL3Transform::Sector2Slice(slice,padrow,cursec,currow); |
e1a32fa5 | 985 | if(slice<fMinSlice || slice>fMaxSlice) continue; |
986 | AliL3SpacePointData *points = fClusters[slice][0]; | |
1b76fd18 | 987 | |
988 | Int_t index = fRowid[slice][padrow]; | |
989 | if(!fDigitsTree->GetEvent(index)) | |
990 | printf("AliL3Evaluate::EvaluatePoints : ERROR IN DIGITSTREE\n"); | |
294ee877 | 991 | printf("Checking slice %d padrow %d with %d clusters\n",slice,padrow,num_of_offline); |
1b76fd18 | 992 | |
e1a32fa5 | 993 | for(UInt_t c=0; c<fNcl[slice][0]; c++) |
78127e35 | 994 | { |
e1a32fa5 | 995 | if(points[c].fPadRow!=padrow) continue; |
e1a32fa5 | 996 | for(Int_t m=0; m<num_of_offline; m++) |
78127e35 | 997 | { |
e1a32fa5 | 998 | AliComplexCluster *cluster = (AliComplexCluster *)clusters->UncheckedAt(m); |
999 | Int_t mcId = cluster->fTracks[0]; | |
294ee877 | 1000 | //Int_t mcId = cluster->GetLabel(0); |
1001 | if(mcId <0) continue; | |
e1a32fa5 | 1002 | TParticle *part = gAlice->Particle(mcId); |
294ee877 | 1003 | |
e1a32fa5 | 1004 | Float_t xyz_cl[3] = {points[c].fX,points[c].fY,points[c].fZ}; |
294ee877 | 1005 | Float_t xyz_ex[3]; |
494fad94 | 1006 | AliL3Transform::Global2Raw(xyz_cl,cursec,currow); |
e170e0f9 | 1007 | if(fDigits->GetTrackID((Int_t)rint(xyz_cl[2]),(Int_t)rint(xyz_cl[1]),0)!=mcId && |
1008 | fDigits->GetTrackID((Int_t)rint(xyz_cl[2]),(Int_t)rint(xyz_cl[1]),1)!=mcId && | |
1009 | fDigits->GetTrackID((Int_t)rint(xyz_cl[2]),(Int_t)rint(xyz_cl[1]),2)!=mcId) | |
e1a32fa5 | 1010 | continue; |
494fad94 | 1011 | AliL3Transform::Raw2Local(xyz_ex,cursec,currow,cluster->fY,cluster->fX); |
1012 | AliL3Transform::Raw2Local(xyz_cl,cursec,currow,xyz_cl[1],xyz_cl[2]); | |
294ee877 | 1013 | Float_t resy = xyz_cl[1] - xyz_ex[1];//cluster->GetY() |
1014 | Float_t resz = xyz_cl[2] - xyz_ex[2];//cluster->GetZ() | |
e1a32fa5 | 1015 | |
294ee877 | 1016 | ntuppel->Fill(slice,padrow,resy,resz,xyz_ex[2],part->Pt()); |
78127e35 | 1017 | } |
e1a32fa5 | 1018 | } |
294ee877 | 1019 | arr->ClearRow(cursec,currow); |
78127e35 | 1020 | } |
e1a32fa5 | 1021 | |
78127e35 | 1022 | return ntuppel; |
1023 | } | |
355debe1 | 1024 | |
1025 | void AliL3Evaluate::GetCFeff(Char_t *outfile) | |
78127e35 | 1026 | { |
78127e35 | 1027 | |
355debe1 | 1028 | TNtuple *ntuppel = new TNtuple("ntuppel","Cluster finder efficiency","row:ncrossings:nclusters"); |
1029 | ||
1030 | AliTPC *TPC = (AliTPC*)gAlice->GetDetector("TPC"); | |
1031 | ||
1032 | TPC->SetParam(fParam); | |
78127e35 | 1033 | |
355debe1 | 1034 | Int_t ver = TPC->IsVersion(); |
494fad94 | 1035 | LOG(AliL3Log::kInformational,"AliL3Evaluate::GetCFeff","TPC version") |
355debe1 | 1036 | <<"TPC version "<<ver<<" found on file"<<ENDLOG; |
1037 | ||
1038 | Int_t zero=TPC->GetParam()->GetZeroSup(); | |
1039 | ||
1040 | Int_t np = gAlice->GetNtrack(); | |
1041 | ||
1042 | Int_t crossed,recs; | |
1043 | Int_t *count = new Int_t[np]; //np number of particles. | |
1044 | Int_t i; | |
1045 | Float_t xyz[3]; | |
1046 | for (i=0; i<np; i++) count[i]=0; | |
1047 | for(Int_t sl=fMinSlice; sl<=fMaxSlice; sl++) | |
78127e35 | 1048 | { |
355debe1 | 1049 | for (i=0; i<=175; i++) |
78127e35 | 1050 | { |
355debe1 | 1051 | crossed=0; |
1052 | recs=0; | |
1053 | Int_t index = fRowid[sl][i]; | |
1054 | if (!fDigitsTree->GetEvent(index)) continue; | |
1055 | Int_t sec,row; | |
1056 | fParam->AdjustSectorRow(fDigits->GetID(),sec,row); | |
1057 | fDigits->First(); | |
1058 | do { | |
1059 | Int_t it=fDigits->CurrentRow(), ip=fDigits->CurrentColumn(); | |
1060 | Short_t dig = fDigits->GetDigit(it,ip); | |
1061 | ||
1062 | if(dig<=fParam->GetZeroSup()) continue; | |
1063 | if(it < fParam->GetMaxTBin()-1 && it > 0) | |
1064 | if(fDigits->GetDigit(it+1,ip) <= fParam->GetZeroSup() | |
1065 | && fDigits->GetDigit(it-1,ip) <= fParam->GetZeroSup()) | |
1066 | continue; | |
1067 | ||
494fad94 | 1068 | AliL3Transform::Raw2Local(xyz,sec,row,ip,it); |
355debe1 | 1069 | if(fParam->GetPadRowRadii(sec,row)<230./250.*fabs(xyz[2])) |
1070 | continue; | |
1071 | ||
1072 | ||
1073 | Int_t idx0=fDigits->GetTrackID(it,ip,0); | |
1074 | Int_t idx1=fDigits->GetTrackID(it,ip,1); | |
1075 | Int_t idx2=fDigits->GetTrackID(it,ip,2); | |
1076 | ||
1077 | if (idx0>=0 && dig>=zero) count[idx0]+=1; | |
1078 | if (idx1>=0 && dig>=zero) count[idx1]+=1; | |
1079 | if (idx2>=0 && dig>=zero) count[idx2]+=1; | |
1080 | } while (fDigits->Next()); | |
1081 | ||
1082 | for (Int_t j=0; j<np; j++) | |
78127e35 | 1083 | { |
355debe1 | 1084 | if (count[j]>1) //at least two digits at this padrow |
494fad94 | 1085 | { |
1086 | crossed++; | |
1087 | count[j]=0; | |
1088 | } | |
355debe1 | 1089 | } |
1090 | AliL3SpacePointData *points = fClusters[sl][0]; | |
494fad94 | 1091 | for(UInt_t k=0; k<fNcl[sl][0]; k++) |
355debe1 | 1092 | { |
494fad94 | 1093 | if(points[k].fPadRow!=i) continue; |
355debe1 | 1094 | recs++; |
78127e35 | 1095 | } |
355debe1 | 1096 | ntuppel->Fill(i,crossed,recs); |
78127e35 | 1097 | } |
494fad94 | 1098 | |
78127e35 | 1099 | } |
355debe1 | 1100 | delete[] count; |
1101 | ||
1102 | TFile *file = TFile::Open(outfile,"RECREATE"); | |
1103 | ntuppel->Write(); | |
1104 | file->Close(); | |
1105 | ||
78127e35 | 1106 | } |
355debe1 | 1107 | |
78127e35 | 1108 | Bool_t AliL3Evaluate::GetParticleCrossingPoint(TParticle *part,Int_t slice,Int_t padrow,Float_t *xyz) |
1109 | { | |
1110 | //Calcluate the crossing point between a generated particle and given padrow. | |
1111 | ||
eeddc64d | 1112 | Double_t kappa = AliL3Transform::GetBField()*0.0029980/part->Pt(); |
78127e35 | 1113 | |
1114 | Double_t radius = 1/fabs(kappa); | |
1115 | if(part->GetPdgCode() > 0) kappa = -kappa; | |
1116 | ||
23908b9b | 1117 | Int_t charg; |
1118 | if(kappa>0) | |
1119 | charg=-1; | |
1120 | else | |
1121 | charg=1; | |
1122 | ||
78127e35 | 1123 | Float_t angl[1] = {part->Phi()}; |
1124 | ||
494fad94 | 1125 | AliL3Transform::Global2LocalAngle(angl,slice); |
78127e35 | 1126 | |
1127 | Double_t charge = -1.*kappa; | |
84122e15 | 1128 | Double_t trackPhi0 = angl[0] + charge*0.5*AliL3Transform::Pi()/fabs(charge); |
78127e35 | 1129 | |
1130 | Double_t x0=0; | |
1131 | Double_t y0=0; | |
1132 | Double_t xc = x0 - radius * cos(trackPhi0); | |
1133 | Double_t yc = y0 - radius * sin(trackPhi0); | |
1134 | ||
1135 | //printf("radius %f xc %f yc %f\n",radius,xc,yc); | |
1136 | ||
494fad94 | 1137 | Double_t xHit = AliL3Transform::Row2X(padrow); |
78127e35 | 1138 | xyz[0] = xHit; |
1139 | Double_t aa = (xHit - xc)*(xHit - xc); | |
1140 | Double_t r2 = radius*radius; | |
1141 | if(aa > r2) | |
1142 | return false; | |
1143 | ||
1144 | Double_t aa2 = sqrt(r2 - aa); | |
1145 | Double_t y1 = yc + aa2; | |
1146 | Double_t y2 = yc - aa2; | |
1147 | xyz[1] = y1; | |
1148 | if(fabs(y2) < fabs(y1)) xyz[1] = y2; | |
1149 | ||
23908b9b | 1150 | Double_t tgl = part->Pz()/part->Pt(); |
1151 | ||
1152 | Double_t yHit = xyz[1]; | |
1153 | Double_t angle1 = atan2((yHit - yc),(xHit - xc)); | |
84122e15 | 1154 | if(angle1 < 0) angle1 += 2.*AliL3Transform::Pi(); |
23908b9b | 1155 | Double_t angle2 = atan2((0 - yc),(0 - xc)); |
84122e15 | 1156 | if(angle2 < 0) angle2 += 2.*AliL3Transform::Pi(); |
23908b9b | 1157 | Double_t diff_angle = angle1 - angle2; |
84122e15 | 1158 | diff_angle = fmod(diff_angle,2*AliL3Transform::Pi()); |
1159 | if((charg*diff_angle) > 0) diff_angle = diff_angle - charg*2.*AliL3Transform::Pi(); | |
23908b9b | 1160 | Double_t s_tot = fabs(diff_angle)*radius; |
1161 | Double_t zHit = 0 + s_tot*tgl; | |
1162 | xyz[2] = zHit; | |
1163 | ||
78127e35 | 1164 | return true; |
1165 | } | |
1166 |