]>
Commit | Line | Data |
---|---|---|
212b7ec8 | 1 | #include "TObjArray.h" |
2 | #include "TMath.h" | |
3 | ||
4 | #include "AliLog.h" | |
5 | #include "AliMathBase.h" | |
6 | ||
7 | #include "AliTRDseedV1.h" | |
8 | #include "AliTRDcluster.h" | |
9 | #include "AliTRDpadPlane.h" | |
10 | #include "AliTRDtrackletOflHelper.h" | |
11 | ||
12 | ClassImp(AliTRDtrackletOflHelper) | |
13 | //___________________________________________________________________ | |
14 | AliTRDtrackletOflHelper::AliTRDtrackletOflHelper() | |
15 | :TObject() | |
16 | ,fRow(-1) | |
17 | ,fClusters(NULL) | |
18 | ,fPadPlane(NULL) | |
19 | { | |
20 | // Default constructor | |
21 | fCol[0]=144; fCol[1]=0; | |
22 | fTBrange[0]=1; fTBrange[1]=21; | |
23 | } | |
24 | ||
25 | //___________________________________________________________________ | |
26 | AliTRDtrackletOflHelper::AliTRDtrackletOflHelper(const AliTRDtrackletOflHelper &ref) | |
27 | :TObject(ref) | |
28 | ,fRow(-1) | |
29 | ,fClusters(NULL) | |
30 | ,fPadPlane(ref.fPadPlane) | |
31 | { | |
32 | // Copy constructor | |
33 | fRow = ref.fRow; fCol[0]=ref.fCol[0]; fCol[1]=ref.fCol[1]; | |
34 | fTBrange[0] = ref.fTBrange[0];fTBrange[1] = ref.fTBrange[1]; | |
35 | Int_t n(0); | |
36 | if(ref.fClusters){ | |
37 | if((n = ref.fClusters->GetEntriesFast())) { | |
38 | fClusters = new TObjArray(n); | |
39 | for(Int_t ic(n); ic--;) fClusters->AddAt(ref.fClusters->At(ic), ic); | |
40 | } | |
41 | } | |
42 | } | |
43 | ||
44 | //___________________________________________________________________ | |
45 | AliTRDtrackletOflHelper& AliTRDtrackletOflHelper::operator=(const AliTRDtrackletOflHelper &rhs) | |
46 | { | |
47 | if(this != &rhs){ | |
48 | TObject::operator=(rhs); | |
49 | fPadPlane = rhs.fPadPlane; | |
50 | fRow = rhs.fRow; | |
51 | fCol[0] = rhs.fCol[0]; | |
52 | fCol[1] = rhs.fCol[1]; | |
53 | fTBrange[0]= rhs.fTBrange[0]; | |
54 | fTBrange[1]= rhs.fTBrange[1]; | |
55 | if(rhs.fClusters){ | |
56 | Int_t n(rhs.fClusters->GetEntriesFast()); | |
57 | if(!fClusters) fClusters = new TObjArray(n); | |
58 | if(fClusters->GetEntriesFast() != n){ | |
59 | fClusters->Clear(); | |
60 | fClusters->Expand(n); | |
61 | } | |
62 | for(Int_t ic(n); ic--;) fClusters->AddAt(rhs.fClusters->At(ic), ic); | |
63 | } else { | |
64 | if(fClusters) delete fClusters; | |
65 | } | |
66 | } | |
67 | return *this; | |
68 | } | |
69 | ||
70 | //___________________________________________________________________ | |
71 | AliTRDtrackletOflHelper::~AliTRDtrackletOflHelper() | |
72 | { | |
73 | // Clean helper | |
74 | if(fClusters) fClusters->Clear(); | |
75 | delete fClusters; | |
76 | } | |
77 | ||
78 | //___________________________________________________________________ | |
79 | Int_t AliTRDtrackletOflHelper::Expand(TObjArray *cls, Int_t *mark, Int_t groupId) | |
80 | { | |
81 | // Allocate new clusters for this helper. If a subset of clusters is to be allocated | |
82 | // this can be specified via the identifier "groupId" which have to be compared | |
83 | // with individual elements in the array "mark". | |
84 | // | |
85 | // Return total no. of clusters in helper | |
86 | ||
87 | if(!fPadPlane || !fClusters ){ | |
88 | AliError("Helper not initialized."); | |
89 | return 0; | |
90 | } | |
91 | Int_t ncl(fClusters->GetEntriesFast()); | |
92 | Int_t mcl(cls->GetEntriesFast()); | |
93 | fClusters->Expand(mcl + ncl); | |
94 | fCol[0]=144; fCol[1]=0; | |
95 | AliTRDcluster *c(NULL); | |
96 | for(Int_t icl(0), jcl(ncl); icl<mcl; icl++){ | |
97 | if(!(c=(AliTRDcluster*)cls->At(icl))) continue; | |
98 | if(mark && mark[icl]!=groupId) continue; | |
99 | fClusters->AddAt(c, jcl++); | |
100 | Int_t col(c->GetPadCol()); | |
101 | Short_t *sig(c->GetSignals()); | |
102 | for(Int_t icol(0), jcol(col-3); icol<7; icol++, jcol++){ | |
103 | if(sig[icol]==0) continue; | |
104 | if(jcol<fCol[0]) fCol[0]=jcol; | |
105 | if(jcol>fCol[1]) fCol[1]=jcol; | |
106 | } | |
107 | } | |
108 | ||
109 | AliDebug(1, Form("Segment[%d] Clusters[%2d] pad[%3d %3d].", groupId, fClusters->GetEntriesFast(), fCol[0], fCol[1])); | |
110 | return fClusters->GetEntriesFast(); | |
111 | } | |
112 | ||
113 | ||
114 | //___________________________________________________________________ | |
115 | Int_t AliTRDtrackletOflHelper::Init(AliTRDpadPlane *p, TObjArray *cls, Int_t *mark, Int_t groupId) | |
116 | { | |
117 | // Allocate clusters for this helper. If a subset of clusters is to be allocated | |
118 | // this can be specified via the identifier "groupId" which have to be compared | |
119 | // with individual elements in the array "mark". | |
120 | // | |
121 | // Return no. of clusters allocated | |
122 | ||
123 | if(!p){ | |
124 | AliError("PadPlane not initialized."); return 0; | |
125 | } | |
126 | if(!cls){ | |
127 | AliError("Cluster array not initialized."); return 0; | |
128 | } | |
129 | Int_t ncl(cls->GetEntriesFast()); | |
130 | if(ncl<2){ | |
131 | AliDebug(1, Form("Segment[%d] failed n[%d].", groupId, ncl)); return 0; | |
132 | } | |
133 | Int_t mcl(ncl); | |
134 | if(mark){ | |
135 | mcl = 0; | |
136 | for(Int_t icl(ncl); icl--;) | |
137 | if(cls->At(icl) && mark[icl]==groupId) mcl++; | |
138 | } | |
139 | if(mcl<2){ | |
140 | AliDebug(1, Form("Segment[%d] failed n[%d] in group.", groupId, mcl)); return 0; | |
141 | } | |
142 | if(!fClusters) fClusters = new TObjArray(mcl); | |
143 | else{ | |
144 | fClusters->Clear(); | |
145 | fClusters->Expand(mcl); | |
146 | } | |
147 | fCol[0]=144; fCol[1]=0; fRow=-1; | |
148 | AliTRDcluster *c(NULL); | |
149 | for(Int_t icl(0), jcl(0); icl<ncl; icl++){ | |
150 | if(!(c=(AliTRDcluster*)cls->At(icl))) continue; | |
151 | if(mark && mark[icl]!=groupId) continue; | |
152 | fClusters->AddAt(c, jcl++); | |
153 | if(fRow<0) fRow = c->GetPadRow(); | |
154 | Int_t col(c->GetPadCol()); | |
155 | Short_t *sig(c->GetSignals()); | |
156 | for(Int_t icol(0), jcol(col-3); icol<7; icol++, jcol++){ | |
157 | if(sig[icol]==0) continue; | |
158 | if(jcol<fCol[0]) fCol[0]=jcol; | |
159 | if(jcol>fCol[1]) fCol[1]=jcol; | |
160 | } | |
161 | } | |
162 | fPadPlane = p; | |
163 | ||
164 | AliDebug(1, Form("Segment[%d] Clusters[%2d] pad[%3d %3d].", groupId, fClusters->GetEntriesFast(), fCol[0], fCol[1])); | |
165 | return fClusters->GetEntriesFast(); | |
166 | } | |
167 | ||
168 | //___________________________________________________________________ | |
169 | Int_t AliTRDtrackletOflHelper::ClassifyTopology() | |
170 | { | |
171 | // Classify topology and return classification code | |
172 | // 0 - normal tracklet | |
173 | // 1 - delta ray candidate | |
174 | // 2 - secondary candidate | |
175 | // 3 - "elephant" candidate | |
176 | // 4 - unknown topology | |
177 | ||
178 | if(!fClusters){ | |
179 | AliError("Helper not initialized. Missing clusters."); | |
180 | return kUnknown; | |
181 | } | |
182 | Int_t ncl(fClusters->GetEntries()); | |
183 | if(!ncl){ | |
184 | AliError("Helper not initialized. No cluster allocated."); | |
185 | return kUnknown; | |
186 | } | |
187 | const Int_t kRange(22); // DEFINE based on vd | |
188 | ||
189 | // compute local occupancy | |
190 | Int_t localOccupancy[AliTRDseedV1::kNtb]; memset(localOccupancy, 0, AliTRDseedV1::kNtb*sizeof(Int_t)); | |
191 | AliTRDcluster *c(NULL); | |
192 | Double_t sy[kNcls]; Int_t mcl(0); | |
193 | for(Int_t icl(ncl), time; icl--;){ | |
194 | c = (AliTRDcluster*)fClusters->At(icl); | |
195 | time = c->GetPadTime(); | |
196 | if(time==0 || time>=kRange){ | |
197 | sy[icl] = -1; // mark clusters outsde drift volume | |
198 | continue; | |
199 | } | |
200 | // protect against wrong error param. | |
201 | if(c->GetSigmaY2() < 1.e-5) sy[icl] = 0.02; | |
202 | else sy[icl] = TMath::Min(TMath::Sqrt(c->GetSigmaY2()), 0.04); | |
203 | localOccupancy[time]++; | |
204 | mcl++; | |
205 | } | |
206 | if(!mcl){ | |
207 | AliWarning("No clusters in the active area."); | |
208 | return kUnknown; | |
209 | } | |
210 | //compute number of 0 bins and high occupancy | |
211 | Int_t goodOccupancy(0), highOccupancy(0), lowOccupancy(0); | |
212 | for(Int_t itb(1); itb<kRange; itb++){ | |
213 | switch(localOccupancy[itb]){ | |
214 | case 0: lowOccupancy++; break; | |
215 | case 1: goodOccupancy++; break; | |
216 | default: highOccupancy++; break; | |
217 | } | |
218 | } | |
219 | AliDebug(2, Form("H[%2d | %5.2f%%] L[%2d | %5.2f%%] N[%2d | %5.2f%% | %5.2f%%]", | |
220 | highOccupancy, 100.*highOccupancy/kRange, | |
221 | lowOccupancy, 100.*lowOccupancy/kRange, | |
222 | goodOccupancy, 100.*goodOccupancy/kRange, 100.*goodOccupancy/mcl)); | |
223 | ||
224 | // filter | |
225 | Double_t dy[kNcls]; | |
226 | if(goodOccupancy==mcl) return kNormal; | |
227 | else if(Double_t(goodOccupancy)/kRange > 0.9){ | |
228 | if(highOccupancy == 0 ) { | |
229 | if(!FitPSR(dy)) return kUnknown; | |
230 | Int_t nin(0); | |
231 | for(Int_t idy(ncl); idy--;){ | |
232 | if(sy[idy]<0.) continue; | |
233 | if(dy[idy] > 3*sy[idy]) continue; | |
234 | nin ++; | |
235 | } | |
236 | if(Double_t(nin)/mcl > 0.8) return kNormal; | |
237 | else return kUnknown; | |
238 | } else return kDeltaRay; | |
239 | } else return kUnknown; | |
240 | } | |
241 | ||
242 | ||
243 | //___________________________________________________________________ | |
244 | void AliTRDtrackletOflHelper::FindSolidCls(Bool_t *mark, Int_t *q) | |
245 | { | |
246 | // Find clusters produced by large fluctuations of energy deposits | |
247 | // Largest charge and well separation from neighbors | |
248 | ||
0a62661e | 249 | Int_t ntb(AliTRDseedV1::kNtb); |
250 | Int_t idx[ntb+1]; | |
251 | TMath::Sort(ntb, q, idx, kTRUE); | |
212b7ec8 | 252 | Int_t qmax = Int_t(0.3*q[idx[0]]); |
253 | mark[0] = kFALSE; | |
0a62661e | 254 | for(Int_t icl(ntb-5); icl<ntb; icl++) mark[icl] = kFALSE; |
255 | for(Int_t icl(0); icl<ntb; icl++){ | |
212b7ec8 | 256 | Int_t jcl(idx[icl]); |
257 | if(!mark[jcl]) continue; | |
258 | if(q[jcl-1]>q[jcl] || q[jcl+1]>q[jcl]){ | |
259 | mark[jcl] = kFALSE; | |
260 | continue; | |
261 | } | |
262 | if(q[jcl] < qmax){ | |
263 | mark[jcl] = kFALSE; | |
264 | continue; | |
265 | } | |
266 | for(Int_t kcl=TMath::Max(0, jcl-2); kcl<jcl+3; kcl++){ | |
267 | if(kcl==jcl) continue; | |
268 | mark[kcl] = kFALSE; | |
269 | } | |
270 | } | |
271 | } | |
272 | ||
273 | //___________________________________________________________________ | |
274 | Bool_t AliTRDtrackletOflHelper::FitPSR(Double_t dy[200], Bool_t useSolid) | |
275 | { | |
276 | // Fit tracklet in Pad System of Reference [PSR] to avoid uncertainty related to | |
277 | // Lorentz angle correction | |
278 | ||
279 | Bool_t mark[200]; memset(mark, 0, 200*sizeof(Bool_t)); | |
280 | Int_t q[200]; memset(q, 0, 200*sizeof(Int_t)); | |
281 | Int_t ncl(fClusters->GetEntries()); | |
282 | AliTRDcluster *c(NULL); | |
283 | for(Int_t icl(ncl); icl--;){ | |
284 | c = (AliTRDcluster*)fClusters->At(icl); | |
285 | if(c->GetPadRow() != fRow) continue; | |
286 | Int_t time(c->GetPadTime()); | |
287 | mark[time] = kTRUE; q[time] = Int_t(c->GetQ()); | |
288 | } | |
289 | if(useSolid) FindSolidCls(mark, q); | |
290 | ||
291 | Double_t | |
292 | x[AliTRDseedV1::kNtb], y[AliTRDseedV1::kNtb], sy[AliTRDseedV1::kNtb], | |
293 | xf[AliTRDseedV1::kNtb], yf[AliTRDseedV1::kNtb], syf[AliTRDseedV1::kNtb], | |
294 | par[3]; | |
295 | Int_t jcl(0), kcl(0); | |
296 | for(Int_t icl(0); icl<ncl; icl++){ | |
297 | c = (AliTRDcluster*)fClusters->At(icl); | |
298 | if(c->GetPadRow() != fRow) continue; | |
299 | Int_t col(c->GetPadCol()); | |
300 | Int_t time(c->GetPadTime()); | |
301 | Double_t center(c->GetCenter()); | |
302 | Double_t cw(fPadPlane->GetColSize(col)); | |
303 | //Double_t corr = AliTRDcluster::GetYcorr(AliTRDgeometry::GetLayer(det), center); | |
304 | y[jcl] = fPadPlane->GetColPos(col) + (.5 + center)*cw /*+ corr*/; | |
305 | x[jcl] = c->GetX(); | |
306 | sy[jcl]= TMath::Sqrt(c->GetSigmaY2()); | |
307 | if(mark[time]){ | |
308 | yf[kcl] = y[jcl]; | |
309 | xf[kcl] = x[jcl]; | |
310 | syf[kcl]= useSolid ? 0.5/time:sy[jcl]; | |
311 | kcl++; | |
312 | } | |
313 | jcl++; | |
314 | } | |
315 | Fit(kcl, xf, yf, syf, par); | |
316 | ||
317 | for(Int_t icl(0); icl<jcl; icl++){ | |
318 | Double_t dx(x[icl] - par[2]); | |
319 | dy[icl] = y[icl] - (par[0] + par[1]*dx); | |
320 | } | |
321 | return kTRUE; | |
322 | } | |
323 | ||
324 | //___________________________________________________________________ | |
325 | Bool_t AliTRDtrackletOflHelper::Fit(Int_t n, Double_t *x, Double_t *y, Double_t *sy, Double_t *par, Double_t sCut, Double_t *cov) | |
326 | { | |
327 | // Iterative robust tracklet fit | |
328 | if(n<3) return kFALSE; | |
2f4384e6 | 329 | //select reference radial position |
330 | if(par[2]<0.){ //compute reference radial position as <x> | |
212b7ec8 | 331 | par[2] = 0.; for(Int_t ic(n); ic--;) par[2] += x[ic]; par[2] /= n; |
332 | } | |
333 | AliTRDtrackerV1::AliTRDLeastSquare &f=Fitter(); | |
334 | for(Int_t iter(0); iter<3; iter++){ | |
335 | f.Reset(); | |
336 | Int_t jp(0); | |
337 | for(Int_t ip(0); ip<n; ip++){ | |
338 | Double_t dx(x[ip]-par[2]); | |
339 | Double_t dy(y[ip] - (par[0] + par[1]*dx)); | |
340 | if(iter && TMath::Abs(dy)>sCut*sy[ip]) continue; | |
341 | f.AddPoint(&dx, y[ip], sy[ip]); | |
342 | jp++; | |
343 | } | |
344 | if(jp<3) continue; | |
345 | if(!f.Eval()) continue; | |
346 | par[0]=f.GetFunctionParameter(0); | |
347 | par[1]=f.GetFunctionParameter(1); | |
348 | if(cov) f.GetCovarianceMatrix(cov); | |
349 | AliDebugGeneral("AliTRDtrackletOflHelper::Fit()", 2, Form("Iter[%d] Ncl[%2d/%2d] par[%f %f %f]", iter, jp, n, par[0], par[1], par[2])); | |
350 | } | |
351 | return kTRUE; | |
352 | } | |
353 | ||
354 | //___________________________________________________________________ | |
355 | Bool_t AliTRDtrackletOflHelper::Fit(Double_t *par, Double_t sCut) const | |
356 | { | |
357 | // Wrapper for clusters attach to this for static Fit function | |
358 | // | |
359 | Int_t n(fClusters->GetEntriesFast()), jc(0), dr(0); | |
360 | Double_t corr = TMath::Tan(TMath::DegToRad()*fPadPlane->GetTiltingAngle())* | |
361 | fPadPlane->GetLengthIPad(); | |
362 | Double_t x[kNcls], y[kNcls], sy[kNcls]; | |
363 | AliTRDcluster *c(NULL); | |
364 | for(Int_t ic(0); ic<n; ic++){ | |
365 | c = (AliTRDcluster*)fClusters->At(ic); | |
366 | if(!c->IsInChamber()) continue; | |
367 | dr = c->GetPadRow() - fRow; | |
368 | x[jc] = c->GetX(); | |
369 | y[jc] = c->GetY()+corr*dr; | |
370 | sy[jc]= c->GetSigmaY2()>0?(TMath::Min(TMath::Sqrt(c->GetSigmaY2()), 0.08)):0.08; | |
371 | jc++; | |
372 | } | |
373 | return Fit(jc, x, y, sy, par, sCut); | |
374 | } | |
375 | ||
376 | //___________________________________________________________________ | |
377 | void AliTRDtrackletOflHelper::GetColSignals(Int_t col, Int_t adc[32], Bool_t mainRow) const | |
378 | { | |
379 | memset(adc, 0, 32*sizeof(Int_t)); | |
380 | if(col<fCol[0] || col>fCol[1]) return; | |
381 | ||
382 | AliTRDcluster *cc(NULL); | |
383 | for(Int_t ic(fClusters->GetEntriesFast()); ic--;){ | |
384 | cc = (AliTRDcluster*)fClusters->At(ic); | |
385 | if((mainRow && cc->GetPadRow()!=fRow) || | |
386 | (!mainRow && cc->GetPadRow()==fRow)) continue; | |
387 | Short_t *sig = cc->GetSignals(); | |
388 | Int_t padcol(cc->GetPadCol()); | |
389 | Int_t time(cc->GetPadTime()); | |
390 | for(Int_t icol(0), jcol(padcol-3); icol<7; icol++, jcol++){ | |
391 | if(jcol!=col) continue; | |
392 | adc[time]+=sig[icol]; | |
393 | } | |
394 | } | |
395 | } | |
396 | ||
397 | //___________________________________________________________________ | |
2f4384e6 | 398 | Int_t AliTRDtrackletOflHelper::GetRMS(Double_t &r, Double_t &m, Double_t &s, Double_t xm) const |
212b7ec8 | 399 | { |
2f4384e6 | 400 | // Calculate Rotation[r], Mean y[m] (at radial position [xm]) and Sigma[s] (of a gaussian distribution in the tracklet [SR]) |
212b7ec8 | 401 | // for clusters attach to this helper. The Rotation and Mean are calculated without tilt correction option. |
402 | // It returns the number of clusters in 1 sigma cut. | |
403 | ||
404 | Int_t n(fClusters->GetEntriesFast()); | |
2f4384e6 | 405 | if(n<=2) return n; |
406 | Double_t par[3] = {0., 0., -1.}; | |
407 | if(xm>0.) par[2] = xm; // select reference radial position from outside | |
212b7ec8 | 408 | Fit(par); |
2f4384e6 | 409 | m = par[0]; |
212b7ec8 | 410 | r = par[1]; |
2f4384e6 | 411 | xm= par[2]; |
212b7ec8 | 412 | |
413 | Double_t corr = TMath::Tan(TMath::DegToRad()*fPadPlane->GetTiltingAngle())* | |
414 | fPadPlane->GetLengthIPad(); | |
415 | Double_t y[kNcls]; | |
416 | AliTRDcluster *c(NULL); | |
417 | for(Int_t ic(n); ic--;){ | |
418 | c = (AliTRDcluster*)fClusters->At(ic); | |
419 | Double_t x(c->GetX() - xm); | |
420 | Int_t dr(c->GetPadRow() - fRow); | |
421 | y[ic] = c->GetY()+corr*dr - (par[0] + par[1]*x); | |
422 | } | |
423 | Double_t m1(0.); | |
424 | AliMathBase::EvaluateUni(n, y, m1, s, 0); | |
212b7ec8 | 425 | Int_t n0(0); |
426 | for(Int_t ic(n); ic--;){ | |
427 | c = (AliTRDcluster*)fClusters->At(ic); | |
428 | Double_t sy = c->GetSigmaY2()>0?(TMath::Min(TMath::Sqrt(c->GetSigmaY2()), 0.08)):0.08; | |
429 | if(TMath::Abs(y[ic]-m1) <= sy) n0++; | |
430 | } | |
431 | return n0; | |
432 | } | |
433 | ||
434 | //___________________________________________________________________ | |
2f4384e6 | 435 | Double_t AliTRDtrackletOflHelper::GetQ() const |
436 | { | |
437 | // Calculate total charge NOT normalized to inclination | |
438 | ||
439 | Double_t q(0.); | |
440 | Int_t n(fClusters->GetEntriesFast()); | |
441 | AliTRDcluster *c(NULL); | |
442 | for(Int_t ic(n); ic--;){ | |
443 | c = (AliTRDcluster*)fClusters->At(ic); | |
444 | q += TMath::Abs(c->GetQ()); | |
445 | } | |
446 | return q; | |
447 | } | |
448 | ||
449 | //___________________________________________________________________ | |
212b7ec8 | 450 | Double_t AliTRDtrackletOflHelper::GetSyMean() const |
451 | { | |
452 | // Calculate mean uncertainty of clusters | |
453 | ||
454 | Double_t sym(0.); | |
455 | Int_t n(fClusters->GetEntriesFast()); | |
456 | AliTRDcluster *c(NULL); | |
457 | for(Int_t ic(n); ic--;){ | |
458 | c = (AliTRDcluster*)fClusters->At(ic); | |
459 | sym += c->GetSigmaY2()>0?(TMath::Min(TMath::Sqrt(c->GetSigmaY2()), 0.08)):0.08; | |
460 | } | |
461 | return sym/=n; | |
462 | } | |
463 | ||
464 | //___________________________________________________________________ | |
465 | Int_t AliTRDtrackletOflHelper::Segmentation(Int_t n, Double_t *x, Double_t *y, Int_t *Index) | |
466 | { | |
467 | // Segmentation of clusters in the tracklet roads | |
468 | // | |
469 | // The user supply the coordinates of the clusters (x,y) and their number "n". Also | |
470 | // The array "Index" has to be allocated by the user with a size equal or larger than "n". | |
471 | // On return the function returns the number of segments found and the array "Index" i-th element | |
472 | // is filled with the index of the tracklet segment to which the i-th cluster was assigned too. | |
473 | // | |
474 | // Observation | |
475 | // The parameter which controls the segmentation is set inside the function "kGapSize" and it is used | |
476 | // for both "x" and "y" segmentations. An improvement can be a parametrization as function of angle of | |
477 | // incidence. | |
478 | // | |
479 | // author: | |
480 | // Alex Bercuci <a.bercuci@gsi.de> | |
481 | ||
482 | if(!n || !x || !y || !Index){ | |
483 | AliErrorGeneral("AliTRDtrackletOflHelper::Segmentation()", "One of the input arrays non initialized."); | |
484 | return 0; | |
485 | } | |
e08c6080 | 486 | const Int_t kBuffer = 200; |
487 | if(n>kBuffer){ | |
488 | AliWarningGeneral("AliTRDtrackletOflHelper::Segmentation()", Form("Input array size %d exceed buffer %d. Truncate.", n, kBuffer)); | |
489 | n = kBuffer; | |
490 | } | |
212b7ec8 | 491 | const Double_t kGapSize(0.2); // cm |
492 | Int_t ng(0), | |
493 | nc(0); | |
e08c6080 | 494 | Double_t xx[kBuffer], dy; |
0a62661e | 495 | Int_t idx[kBuffer+1], jdx[kBuffer], kdx[kBuffer]; |
212b7ec8 | 496 | TMath::Sort(n, y, idx); |
497 | for(Int_t iy(0); iy<n; iy++){ | |
498 | dy = iy>0?(TMath::Abs(y[idx[iy-1]]-y[idx[iy]])):0.; | |
499 | if(dy>kGapSize){ | |
500 | TMath::Sort(nc, xx, jdx); | |
501 | for(Int_t ic(0), jc0, jc1; ic<nc; ic++){ | |
502 | jc0 = ic>0?kdx[jdx[ic-1]]:0; | |
503 | jc1 = kdx[jdx[ic]]; | |
504 | dy = TMath::Abs(y[jc0] - y[jc1]); | |
505 | if(ic && dy>kGapSize) ng++; | |
506 | Index[jc1] = ng; | |
507 | AliDebugGeneral("AliTRDtrackletOflHelper::Segmentation()", 4, Form(" y[%2d]=%+f x[%+f] %2d -> %2d -> %2d ng[%2d]", jc1, y[jc1], x[jc1], ic, jdx[ic], jc1, ng)); | |
508 | } | |
509 | ng++; | |
510 | nc=0; | |
511 | } | |
512 | xx[nc] = x[idx[iy]]; | |
513 | kdx[nc]= idx[iy]; | |
514 | AliDebugGeneral("AliTRDtrackletOflHelper::Segmentation()", 4, Form("y[%2d]=%+f -> %2d", idx[iy], y[idx[iy]], nc)); | |
515 | nc++; | |
516 | } | |
517 | if(nc){ | |
518 | TMath::Sort(nc, xx, jdx); | |
519 | for(Int_t ic(0), jc0, jc1; ic<nc; ic++){ | |
520 | jc0 = ic>0?kdx[jdx[ic-1]]:0; | |
521 | jc1 = kdx[jdx[ic]]; | |
522 | dy = TMath::Abs(y[jc0] - y[jc1]); | |
523 | if(ic && dy>kGapSize) ng++; | |
524 | Index[jc1] = ng; | |
525 | AliDebugGeneral("AliTRDtrackletOflHelper::Segmentation()", 4, Form(" y[%2d]=%+f x[%+f|%+f] %2d -> %2d -> %2d ng[%2d]\n", jc1, y[jc1], xx[jdx[ic]], x[jc1], ic, jdx[ic], jc1, ng)); | |
526 | } | |
527 | ng++; | |
528 | } | |
529 | return ng; | |
530 | } | |
531 | ||
532 | //___________________________________________________________________ | |
533 | void AliTRDtrackletOflHelper::SetTbRange(Float_t t0, Float_t vd) | |
534 | { | |
535 | // Set first time bin and total number of time bins corresponding to clusters | |
536 | // in chamber based on the calibrated info "t0" and drift velocity "vd" | |
537 | ||
538 | // TO CHECK | |
539 | fTBrange[0] = Int_t(t0*0.1); | |
540 | fTBrange[1] = Int_t(vd*0.1); | |
541 | } | |
542 | ||
543 | #include "TH2.h" | |
544 | #include "TGraph.h" | |
545 | #include "TGraphErrors.h" | |
546 | #include "TCanvas.h" | |
547 | #include "TLegend.h" | |
548 | #include "TGaxis.h" | |
549 | //___________________________________________________________________ | |
550 | void AliTRDtrackletOflHelper::View(TVirtualPad *vpad) | |
551 | { | |
552 | // Visualization support. Draw this tracklet segment | |
553 | ||
554 | ||
555 | Int_t row(-1), det(-1); | |
556 | Int_t n(fClusters->GetEntriesFast()); | |
557 | AliInfo(Form("Processing Clusters[%2d] Class[%d].", n ,ClassifyTopology())); | |
558 | ||
559 | // prepare drawing objects | |
560 | Int_t ncols(fPadPlane->GetNcols()); | |
561 | Double_t cw(fPadPlane->GetColSize(1)); | |
562 | TH2 *h2 = new TH2I("h2pm", ";y_{Local} [cm]; pad time;Charge", ncols, fPadPlane->GetColPos(1)-cw, fPadPlane->GetColPos(ncols-1)+cw, 31, -30.5, 0.5); | |
563 | h2->SetMarkerColor(kWhite);h2->SetMarkerSize(2.); | |
564 | h2->SetFillColor(9); | |
565 | ||
566 | TGraph *gcls = new TGraph(n); | |
567 | gcls->SetMarkerColor(kBlack);gcls->SetMarkerStyle(20); | |
568 | TGraph *gclsLC = new TGraph(n); | |
569 | gclsLC->SetMarkerColor(kBlack);gclsLC->SetMarkerStyle(28); | |
570 | TGraph *gclsSC = new TGraph(n); | |
571 | gclsSC->SetMarkerColor(kBlack);gclsSC->SetMarkerStyle(4);gclsSC->SetMarkerSize(1.5); | |
572 | TGraph *gFP = new TGraph(n); | |
573 | gFP->SetMarkerColor(kBlack);gFP->SetLineWidth(2); | |
574 | ||
575 | // fill signal data | |
576 | Bool_t map[200]; memset(map, 0, 200*sizeof(Bool_t)); | |
577 | Int_t qa[200]; memset(qa, 0, 200*sizeof(Int_t)); | |
578 | AliTRDcluster *cc(NULL); Int_t tm = 30; Double_t ym=0.; | |
579 | for(Int_t ic(n); ic--;){ | |
580 | cc = (AliTRDcluster*)fClusters->At(ic); | |
581 | if(cc->GetPadRow() != fRow) continue; //TODO extend for 2 pad rows | |
582 | Short_t *sig = cc->GetSignals(); | |
583 | Int_t col(cc->GetPadCol()); | |
584 | det = cc->GetDetector(); | |
585 | row = cc->GetPadRow(); | |
586 | Int_t time(cc->GetPadTime()); | |
587 | map[time] = kTRUE; qa[time] = (Int_t)cc->GetQ(); | |
588 | for(Int_t ipad(0), jpad(col-2); ipad<7; ipad++, jpad++){ | |
589 | Int_t q = (Int_t)h2->GetBinContent(jpad, 31-time); | |
590 | h2->SetBinContent(jpad, 31-time, q+sig[ipad]); | |
591 | } | |
592 | Double_t y0 = fPadPlane->GetColPos(col) + (.5 + cc->GetCenter())*cw; | |
593 | //h2->GetXaxis()->GetBinCenter(pad)+cc->GetCenter(); | |
594 | gcls->SetPoint(ic, y0, -time); | |
595 | if(time <= tm) {ym = cc->GetY() - y0; tm = time;} | |
596 | } | |
597 | ||
598 | // draw special clusters (solid and Lorentz corrected and fit) | |
599 | FindSolidCls(map, qa); | |
600 | Double_t ddy[200]; FitPSR(ddy, kTRUE); | |
601 | Double_t dt, dy; | |
602 | for(Int_t ic(0), jc(0); ic<n; ic++){ | |
603 | cc = (AliTRDcluster*)fClusters->At(ic); | |
604 | if(cc->GetPadRow() != fRow) continue; //TODO extend for 2 pad rows | |
605 | gcls->GetPoint(ic, dy, dt); | |
606 | gclsLC->SetPoint(ic, cc->GetY()-ym, dt); | |
607 | gFP->SetPoint(ic, dy-ddy[ic], dt); | |
608 | if(map[cc->GetPadTime()]) gclsSC->SetPoint(jc++, dy, dt); | |
609 | } | |
610 | ||
611 | ||
612 | // prepare frame histo | |
613 | h2->SetName(Form("h2s%03d%02d", det, row)); | |
614 | Int_t binSrt(fCol[0]+1), binSop(fCol[1]+1); | |
615 | TH1 *h1(NULL); | |
616 | ||
617 | // show everything | |
618 | if(!vpad){ | |
619 | vpad = new TCanvas(Form("c%03d%02d", det, row), Form("D %03d [%02d_%d_%d] R[%02d]", det, AliTRDgeometry::GetSector(det), AliTRDgeometry::GetStack(det), AliTRDgeometry::GetLayer(det), row), 700, 500); | |
620 | } | |
621 | TVirtualPad *pp(NULL); | |
622 | vpad->Divide(2,1,2.e-5,2.e-5); | |
623 | pp = vpad->cd(1); pp->SetRightMargin(0.0001);pp->SetTopMargin(0.1);pp->SetBorderMode(0); pp->SetFillColor(kWhite); | |
624 | h1 = h2->ProjectionY(); h1->GetYaxis()->SetTitle("Total Charge"); | |
625 | h1->SetTitle(Form("D%03d[%02d_%d_%d] R[%02d]", | |
626 | det,AliTRDgeometry::GetSector(det), AliTRDgeometry::GetStack(det), AliTRDgeometry::GetLayer(det), row)); | |
627 | TAxis *ax(h1->GetXaxis()); for(Int_t ib(0), jb(1); ib<ax->GetNbins(); ib++, jb++) ax->SetBinLabel(jb, Form("%d", -Int_t(ax->GetBinCenter(jb)))); | |
628 | h1->Draw("hbar3"); | |
629 | ||
630 | pp = vpad->cd(2); | |
631 | pp->SetRightMargin(0.15);pp->SetLeftMargin(0.0001);pp->SetTopMargin(0.1); | |
632 | pp->SetBorderMode(0); pp->SetFillColor(kWhite); | |
633 | ax=h2->GetXaxis(); | |
634 | ax->SetRange(TMath::Max(1, binSrt-1), TMath::Min(ncols, binSop+1)); | |
635 | h2->Draw("coltextz"); | |
636 | gPad->Update(); | |
637 | TGaxis *axis = new TGaxis(gPad->GetUxmin(), | |
638 | gPad->GetUymax(), | |
639 | gPad->GetUxmax(), | |
640 | gPad->GetUymax(), | |
641 | TMath::Max(0, binSrt-2)-0.5, TMath::Min(ncols-1, binSop)+0.5, 510,"-L"); | |
642 | ||
643 | axis->SetNdivisions(103+binSop-binSrt); | |
644 | axis->SetTitle("pad col"); | |
645 | axis->Draw(); | |
646 | ||
647 | TLegend *leg = new TLegend(0.01, 0.1, 0.84, 0.21); | |
648 | leg->SetBorderSize(1); leg->SetFillColor(kYellow-9);leg->SetTextSize(0.04); | |
649 | gcls->Draw("p"); leg->AddEntry(gcls, "Cls. in Pad [SR]", "p"); | |
650 | gclsLC->Draw("p"); leg->AddEntry(gclsLC, "Lorentz Corr. Cls.", "p"); | |
651 | gclsSC->Draw("p"); leg->AddEntry(gclsSC, "Solid Cls.", "p"); | |
652 | gFP->Draw("l"); leg->AddEntry(gFP, "Fit in Pad [SR]", "l"); | |
653 | leg->Draw(); | |
654 | } | |
655 | ||
656 | //___________________________________________________________________ | |
657 | AliTRDtrackerV1::AliTRDLeastSquare& AliTRDtrackletOflHelper::Fitter() | |
658 | { | |
659 | static AliTRDtrackerV1::AliTRDLeastSquare f; | |
660 | return f; | |
661 | } |